I’m reading R. Kreps' paper Parameter uncertainty in (log)normal distributions and trying to figure out how the simulations were done. In order to generate Figure 1, Eqn (2.41) was used. So this is what I did in R:

This is the equation they used: $Z_{eff} = v + z\sqrt{\frac{n(1+v^2)}{w}}$ where

Note: This code is for the case when $n=3$

v <- rt(10000000, 1)/sqrt(3-2)
w <- rchisq(10000000,2)
z <- rnorm(n=10000000, m=0, sd=1) 
z_eff_3 <- v + z * sqrt((3*(1+v*v))/w)

However, I got a completely different graph. I would truly appreciate it if someone can take a look and explain to me what I did wrong.

  • 2
    $\begingroup$ n a sample size and Zeff is defined in Eqn (2.41), but where are v, z, and w defined? Are you sure that z is a realisation from a N(1, 1) random variable? $\endgroup$
    – ocram
    Commented Dec 20, 2012 at 4:35

1 Answer 1


Note: In the original post, I noted that (aside from a minor error which has since been corrected) this was actually a problem with plotting the density, not with the simulation. At the time I wasn't sure what was going on, but I found that truncating the data seemed to solve the problem, so I recommended that. I have since figured out the real problem, so I decided to edit this to suggest a better solution.

With $n=3$, the distribution we are simulating from has very heavy tails. Even though the mean is zero and the inter-quartile range is about 4, in 100,000 draws (or more) we get a few on the order of 1,000,000 in absolute value.

By default, the density function estimates the density at 512 equally-spaced points that completely cover the range of the data (plus a little). In this case, we were only interested in plotting the density over a very small portion of that range--for $x \in [0,5]$. As a result, when we plotted that window, we got what looked like a flat horizontal line, which connected the estimates at points on the order of -2000 and 2000. Not very useful.

My original solution, truncating the data, gave decent looking plots, but wasn't really very satisfying. So I played around with it a bit and figured out the problem and how to solve it: the density function has arguments from and to, which give lower and upper bounds between which is estimates the density. By setting from=0, to=5 we get the following plot, without any truncation: enter image description here

This is both easier and more correct (since the density is distorted slightly by the truncation) than my original solution. The complete code I used (adapted slightly from the code provided in the question) is:

v <- rt(100000, 1)/sqrt(3-2)
w <- rchisq(100000,2)
z <- rnorm(n=100000, m=0, sd=1) 
z_eff_3 <- v + z * sqrt((3*(1+v*v))/w)


Hope that helps.

  • $\begingroup$ Many thanks, Jonathan. Can you please share the R code? Thanks again! $\endgroup$
    – user9292
    Commented Dec 20, 2012 at 6:36
  • $\begingroup$ Here it is.. v <- rt(10000000, 1)/sqrt(3-2) w <- rchisq(10000000,2) z <- rnorm(n=10000000, m=0, sd=1) z_eff_3 <- v + z * sqrt((3*(1+v*v))/w) z_eff_3 <- z_eff_3[which((abs(z_eff_3))<10)] plot(density(z_eff_3),xlim=c(0,5),ylim=c(0,1)) $\endgroup$
    – user9292
    Commented Dec 20, 2012 at 9:30
  • $\begingroup$ I tried to produce the CDF in Fig 3 (for n=3), but again, I got a different graph. This is the R code: F_z_eff_3 <- ecdf(z_eff_3) plot(F_z_eff_3, xlim=c(5,10),ylim=c(0.85,0.99)) Can you please take a look and let me know what's wrong. Many thanks! $\endgroup$
    – user9292
    Commented Dec 20, 2012 at 20:38
  • $\begingroup$ If you're still working with the truncated data it's probably because we truncated pretty close in. You should be able to plot the ecdf of the original (untruncated) data, or, if you want to use the same data to plot the ecdf and the smoothed density estimate, truncating further away (at 100) worked for me. $\endgroup$ Commented Dec 20, 2012 at 20:47
  • $\begingroup$ Edited with a better solution that doesn't require truncation. $\endgroup$ Commented Dec 21, 2012 at 0:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.