monte carlo from lognormal distribution?

Hi,

I know how to select from a normal distribution...for example:

Cons = unifrnd(0,800);  

But how does one select from a lognormal distribution? For data I have the data sample size, the mean, standard deviation, 1st, 5th, 10th, 25th, 50th, 75th, 90th, 95th, 99th percentiles for the distribution along with the max value. I've read the sections on lognormal distribution and I can't quite sort out how to work my data into the example presented.

Thank you in advance for any help!

[John D'Errico]

Then you don't actually have a lognormal distribution, in the sense that you don't know the distribution parameters. All you have a quantiles from what you think is a lognormal. The maximum sample value has no real information content at all, since a lognormal is unbounded.

I would probably take your data and use a lognormal fitting tool to estimate the distribution parameters of a lognormal from your data, NOT from the quantiles.

help lognfit

 LOGNFIT Parameter estimates and confidence intervals for lognormal data.
    PARMHAT = LOGNFIT(X) returns a vector of maximum likelihood estimates 
    PARMHAT(1) = MU and PARMHAT(2) = SIGMA of parameters for a lognormal 
    distribution fitting X.  MU and SIGMA are the mean and standard 
    deviation, respectively, of the associated normal distribution.
 
    [PARMHAT,PARMCI] = LOGNFIT(X) returns 95% confidence intervals for the
    parameter estimates.
 
    [PARMHAT,PARMCI] = LOGNFIT(X,ALPHA) returns 100(1-ALPHA) percent
    confidence intervals for the parameter estimates.
 
    [...] = LOGNFIT(X,ALPHA,CENSORING) accepts a boolean vector of the same
    size as X that is 1 for observations that are right-censored and 0 for
    observations that are observed exactly.
 
    [...] = LOGNFIT(X,ALPHA,CENSORING,FREQ) accepts a frequency vector of
    the same size as X.  FREQ typically contains integer frequencies for
    the corresponding elements in X, but may contain any non-integer
    non-negative values.
 
    [...] = LOGNFIT(X,ALPHA,CENSORING,FREQ,OPTIONS) specifies control
    parameters for the iterative algorithm used to compute ML estimates
    when there is censoring.  This argument can be created by a call to
    STATSET.  See STATSET('lognfit') for parameter names and default values.
 
    Pass in [] for ALPHA, CENSORING, or FREQ to use their default values.
 
    With no censoring, SIGMAHAT is the square root of the unbiased estimate
    of the variance of log(X).  With censoring, SIGMAHAT is the maximum
    likelihood estimate.
 
    See also LOGNCDF, LOGNINV, LOGNLIKE, LOGNPDF, LOGNRND, LOGNSTAT, MLE,
             STATSET.

    Documentation for lognfit
       doc lognfit

I might also take those quantile points, then log them, and then convert to z-scores. If the distribution really is lognormal, then the z-scores would follow a nice linear relationship against the associated probabilities. It would at least give your claim that it is lognormal some credence. But we don't see your data or have any real information to help you better.

[Wesser]

Thank you.

My problem is that I don't have access to the actual raw data....the information I provided in my original comment was taken from a published paper. The authors of that paper also presented the attached plot which shows a histogram of the raw data. I'm interested in the red shaded data that presents as lognormal. Does having this information open up more possibilities?

I found this older question that I was hoping would send me down the right path...but I'm not certain it applies exactly...

[John D'Errico]

No. That does not help you at all. A picture of a distribution does not help at all anyway. Again, you don't have the distribution. You have only a set of (probably approximate) quantiles, taken from a distribution that you think MAY be lognormal. And you have not told us what are those quantile values, so I can't really show you how to recover the distribution parameters from the quantiles, probably using a MLE scheme.

[Wesser]

Ok.

As an example, here is one of the sets of data

[Jeff Miller]

The example screen shot includes the sample mean and SD, and you can estimate the parameters of the lognormal directly from those (6.1871 & 0.5064).

But the lognormal doesn't provide a good fit to that dataset. In the lognormal with your mean and SD, the values at the percentiles in your example are the following, which are pretty far from the sample percentile values in your example: 149.76 211.48 254.2 345.7 486.45 684.51 930.88 1118.9 1580.1

ps. Instead of a screen shot, it would be more helpful to post data as text or code so that it can be copied/pasted.

[Wesser]

Thank you Jeff. I'm going to consult with my advisor on how to proceed. Thanks for your time.

[Image Analyst]

That's exactly why I said if you want to mimic that distribution you need to enter in all the values of the bar chart and not try to fit it to some theoretical distribution like a log-normal. You can draw random numbers from the ACTUAL distribution.

[Wesser]

How would I go about doing what you suggest?

My advisor is asking if it is possible to assign a distribution around the geometric mean (50th % of 589 as proxi)?

% 1 5 10 25 50 75 90 95 99 100

37 80 120 282 589 759 916 1022 1239 1466

[Jeff Miller]

@Wesser To do what @Image Analyst suggests, the first step is to take measurements from the histogram of the distribution that you want to use. For each bar in the histogram, determine the value of X at the center of the bar, and the value of Y at the top of the bar. These are the (X,Y) pairs that @Image Analyst refers to in his answer. Note that the scale of the Y values is not important as these can be normalized later, so the heights of the bars can be measured in whatever unit is most convenient (e.g.,. pixels, mm, etc.)

Answer 1

[Image Analyst]

You can get random samples from either the blue curve or the red curve if you know the heights of those bars.

Just use inverse transform sampling.

https://en.wikipedia.org/wiki/Inverse_transform_sampling

Attached is a demo for Rayleigh.

If you need more help, attach the (x,y) data for the red and blue bar charts.

Answer 2

[John D'Errico]

You have not told me the distribution quantiles yet. So I'll make some up.

P = [1, 5, 10, 25, 50, 75, 90, 95, 99]/100;
Q = logninv(P,2,3);

These are just some made up values. But now we can recover the distribution parameters easily enough. A simple test is to look at this plot:

z0 = norminv(P,0,1);
plot(z0,log(Q),'-o')

If that plot is linear, or reasonably close, then this scheme will work easily enough. In this case of course, the plot is exactly linear, since I constructed the data from a known lognoormal distribution.

P1 = polyfit(norminv(P,0,1),log(Q),1);
Mu = P1(2)

Mu = 2.0000

Sigma = P1(1)

Sigma = 3.0000

Now if you want to generate a sample from that distribution, you can use those distribution parameters.

lognsample = lognrnd(Mu,Sigma,[1,5])

lognsample = 1×5

0.0081 6.4422 109.7892 9.9966 3.1125

If the plot I showed above (for your data) is not linear or close to a straight line, then your distribution may not be lognormal.

[Wesser]

Thank you John. I'm going to consult with my advisor on how to proceed. I appricaite your example.