1 / 13

Latin Hypercube Sampling Example

Latin Hypercube Sampling Example. Jake Blanchard Spring 2010. Example. Z=XY X and Y follow exponential distributions  x =1  y =1/2. Step 1. Divide cdfs into even intervals (vertical axis). Step 2. Now sample a number in each section

orinda
Download Presentation

Latin Hypercube Sampling Example

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Latin Hypercube Sampling Example Jake Blanchard Spring 2010 Uncertainty Analysis for Engineers

  2. Example • Z=XY • X and Y follow exponential distributions • x=1 • y=1/2 Uncertainty Analysis for Engineers

  3. Step 1 • Divide cdfs into even intervals (vertical axis) Uncertainty Analysis for Engineers

  4. Step 2 • Now sample a number in each section • For example, pick a random number between 0.8 and 1.0 and use it to get a random value for x • xx=expinv(rx,mux); Uncertainty Analysis for Engineers

  5. Step 3 • Sort the values (shuffle) • nux=xx(randperm(n)); Uncertainty Analysis for Engineers

  6. An Alternative • Instead of using expinv, we can generate the inverse ourselves • Just take the CDF and solve for x • xx=-mux*log(1-rx); Uncertainty Analysis for Engineers

  7. First Script n=10000000; mux=1; muy=2; x=exprnd(mux,n,1); y=exprnd(muy,n,1); mz=mean(x.*y); error=abs(mz-mux*muy)/mux/muy d=linspace(0,1,n+1); rx=unifrnd(d,d+1/n); ry=unifrnd(d,d+1/n); xx=expinv(rx,mux); yy=expinv(ry,muy); nux=xx(randperm(n)); nuy=yy(randperm(n)); mz=mean(nux.*nuy); error=abs(mz-mux*muy)/mux/muy Uncertainty Analysis for Engineers

  8. Alternative n=10000000; mux=1; muy=2; x=exprnd(mux,n,1); y=exprnd(muy,n,1); mz=mean(x.*y); error=abs(mz-mux*muy)/mux/muy d=linspace(0,1,n+1); rx=unifrnd(d,d+1/n); ry=unifrnd(d,d+1/n); xx=-mux*log(1-rx); yy=-muy*log(1-ry); nux=xx(randperm(n)); nuy=yy(randperm(n)); mz=mean(nux.*nuy); error=abs(mz-mux*muy)/mux/muy Uncertainty Analysis for Engineers

  9. Test • Z=XY • X and Y are beta with mean of 1 and 2, respectively • Use simple Monte Carlo • Then use LHS without sorting • Then use LHS with sorting • N=100,000 • For each case, find mean 100 times and then take standard deviation of results Uncertainty Analysis for Engineers

  10. Case 1 n=100000; ntrials=100; mz=zeros(ntrials,1); for i=1:ntrials x=exprnd(mux,n,1); y=exprnd(muy,n,1); mz(i)=mean(x.*y); end std(mz) Uncertainty Analysis for Engineers

  11. Case 2 d=linspace(0,1,n+1); for i=1:ntrials rx=unifrnd(d,d+1/n); rx=rx(1:end-1); ry=unifrnd(d,d+1/n); ry=ry(1:end-1); x=expinv(rx,mux); y=expinv(ry,muy); mz(i)=mean(x.*y); end std(mz) Uncertainty Analysis for Engineers

  12. Case 3 d=linspace(0,1,n+1); for i=1:ntrials rx=unifrnd(d,d+1/n); rx=rx(1:end-1); ry=unifrnd(d,d+1/n); ry=ry(1:end-1); x=expinv(rx,mux); y=expinv(ry,muy); nux=x(randperm(n)); nuy=y(randperm(n)); mz(i)=mean(nux.*nuy); end std(mz) Uncertainty Analysis for Engineers

  13. Results Uncertainty Analysis for Engineers

More Related