CODEWITHSUNDEEP
roverlay
Titi90
Titi90

Reputation: 139

Overlay a t-distribution to my histogram using R?

How can I overlay a t-density to my histogram using R? Here's my function:

simfun <- function(a=56.25102409,b=1.78977412,c=0.08664925,n=18,x1.sd=18.87671,x2.sd=18.87671,e.sd=18.87671) {
   X1 <- rnorm(n, mean=0, sd=x1.sd)
   X2 <- rnorm(n, mean=0, sd=x2.sd) 
   e <-  rnorm(n, mean=0, sd=e.sd)
   Z <- a+b*X1+c*X2+e 
   data.frame(X1,X2,Z)
}

statfun <- function(samples) {
    coef(lm(Z~X1+X2,data=samples))
}

library(plyr)
B=raply(1000,statfun(simfun()))

(hist(B[,2]))

Upvotes: 0

Views: 4163

Answers (2)

Robert
Robert

Reputation: 5152

As November 2019 the way I found to get this plot is with a new plot. Apparently, using dt with curve assumes the x is normalized. The df parameter could be estimated with fitdistr from MASS package.

fit.t.tc <- fitdistr(B[,2], "t", hessian = TRUE)
(param.t=fit.t.tc$estimate)

dh=hist(B[,2], prob=TRUE)
#to get the scaling correct, then do
curve( dt(x, df=param.t["df"]), add=TRUE, col='blue' ) # ??
par(new = TRUE)
ss=seq(range(dh$mids)[1],range(dh$mids)[2],length.out = 1000)
x=((ss-param.t["m"])/param.t["s"])
plot(x,100*dt(x=x,df=param.t["df"]), type="l",
      col="red", lwd=3,xlab="",axes=F,ylab="")

enter image description here

Upvotes: 1

Greg Snow
Greg Snow

Reputation: 49650

Change the last line to:

hist(B[,2], prob=TRUE)

to get the scaling correct, then do

curve( dt(x, df=15), add=TRUE, col='blue' )

changing the df and color to whatever values you want.

Upvotes: 2

Related Questions