/* This SAS/IML Studio program computes the graphics in http://blogs.sas.com/content/iml/2014/06/13/geometry-lognormal-distribution/ Save as Lognormal.sx. Author: Rick Wicklin Date: 06Jun2014 */ proc iml; /* transformation (m, sqrt(v)) --> (mu, sigma) */ start LN2N(p); m = p[,1]; v = p[,2]##2; phi = sqrt(v + m##2); mu = log(m##2 / phi); sigma = sqrt(log(phi##2/m##2)); return( mu || sigma ); finish; /****************************************/ m = 1:30; /* horizontal grid points */ sqrtV = 1:30; /* vertical grid points */ box = expandgrid(m, sqrtV); /* grid */ fbox = LN2N(box); /* image of grid under xform */ declare DataObject dobj; dobj = DataObject.Create("LN", {"mu" "sigma" "m" "sqrtV"}, fbox || box); /******** grid in domain ***********/ declare ScatterPlot domain; domain = ScatterPlot.Create(dobj, "m", "sqrtV"); domain.SetWindowPosition(50,0,50,50); domain.SetAxisLabel(XAXIS, "Mean (m)"); domain.SetAxisLabel(YAXIS, "Standard Deviation (sqrt(v))"); domain.SetTitleText("Lognormal Parameters", true); domain.ShowObs(false); domain.SetGraphAreaBackgroundColor( WHITE ); domain.SetAspectRatio(1); /* draw grid lines */ domain.DrawUseDataCoordinates(); domain.DrawSetRegion( PLOTBACKGROUND ); /* create N evenly spaced points in [min, max] */ start linspace(min, max, N); return( do(min, max, (max-min)/(N-1)) ); finish; x = linspace(min(m), max(m), 101); /* horz lines */ y = linspace(min(sqrtV), max(sqrtV), 101); /* vert lines */ domain.DrawSetPenColor(BLUE); /* draw horizontal lines */ do i = 1 to ncol(sqrtV); y0 = repeat(sqrtV[i], 1, ncol(x)); domain.DrawLine(x, y0); end; domain.DrawSetPenColor(RED); /* draw vertical lines */ do i = 1 to ncol(m); x0 = repeat(m[i], 1, ncol(y)); domain.DrawLine(x0, y); end; /******** image of grid in range ***********/ declare ScatterPlot range; range = ScatterPlot.Create(dobj, "mu", "sigma"); range.SetWindowPosition(50,50,50,50); range.SetAxisLabel(XAXIS, "Mean (mu)"); range.SetAxisLabel(YAXIS, "Standard Deviation (sigma)"); range.SetTitleText("Corresponding Normal Parameters", true); range.ShowObs(false); range.SetGraphAreaBackgroundColor( WHITE ); range.SetAspectRatio(1); range.DrawUseDataCoordinates(); range.DrawSetRegion( PLOTBACKGROUND ); range.DrawSetPenColor(BLUE); /* draw horizontal lines */ do i = 1 to ncol(sqrtV); y0 = repeat(sqrtV[i], 1, ncol(x)); h = LN2N( x` || y0` ); range.DrawLine(h[,1], h[,2]); end; range.DrawSetPenColor(RED); /* draw vertical lines */ do i = 1 to ncol(m); x0 = repeat(m[i], 1, ncol(y)); v = LN2N( x0` || y` ); range.DrawLine(v[,1], v[,2]); end; /***********************************************/ /* Study the inverse transformation near (mu,sigma)=(3,0.5) */ pause; dobj.SelectObsWhere("mu", WHERE_GE, 2.9); dobj.SelectObsWhere("mu", WHERE_LE, 3.1, WF_SUBSETSELECTION); dobj.SelectObsWhere("sigma", WHERE_GE, 0.4, WF_SUBSETSELECTION); dobj.SelectObsWhere("sigma", WHERE_LE, 0.6, WF_SUBSETSELECTION); /* inverse transformation (mu, sigma) --> (m, sqrt(v)) */ start N2LN(p); mu = p[,1]; sigma = p[,2]; sigma2 = sigma##2; m = exp(mu + sigma2/2); v = (exp(sigma2)-1)#exp(2*mu + sigma2); return( m || sqrt(v) ); finish; /* NLPFDD wants function to return a column vector */ start N2LNcol(x); return( T(N2LN(x)) ); finish; /* linearization of inverse transformation at (3, 0.5) */ x0 = {3 0.5}; call nlpfdd(f, J, JpJ, "N2LNcol", x0, {2 . .}); detJ = det(J); print f, J, detJ; pause; /***********************************************/ /* you can make similar pictures for the inverse transformation */ mu = do(2.9, 3.1, 0.01); /* horizontal grid points */ sigma = do(0.4, 0.6, 0.01); /* vertical grid points */ box = expandgrid(mu, sigma); /* grid */ fbox = N2LN(box); /* image of grid under xform */ declare DataObject dobj; dobj = DataObject.Create("LN", {"mu" "sigma" "m" "sqrtV"}, box || fbox); /******** grid in domain ***********/ declare ScatterPlot domain; domain = ScatterPlot.Create(dobj, "mu", "sigma"); domain.SetWindowPosition(0,50,50,50); domain.SetAxisLabel(XAXIS, "Mean (mu)"); domain.SetAxisLabel(YAXIS, "Standard Deviation (sigma)"); domain.SetTitleText("Normal Parameters", true); domain.ShowObs(false); domain.SetGraphAreaBackgroundColor( WHITE ); domain.SetAspectRatio(1); domain.DrawUseDataCoordinates(); x = linspace(min(mu), max(mu), 101); /* horz lines */ y = linspace(min(sigma), max(sigma), 101); /* vert lines */ domain.DrawSetPenColor(BLUE); /* draw horizontal lines */ do i = 1 to ncol(sigma); y0 = repeat(sigma[i], 1, ncol(x)); domain.DrawLine(x, y0); end; domain.DrawSetPenColor(RED); /* draw vertical lines */ do i = 1 to ncol(mu); x0 = repeat(mu[i], 1, ncol(y)); domain.DrawLine(x0, y); end; /******** image of grid in range ***********/ declare ScatterPlot range; range = ScatterPlot.Create(dobj, "m", "sqrtV"); range.SetWindowPosition(0,0,50,50); range.SetAxisLabel(XAXIS, "Mean (m)"); range.SetAxisLabel(YAXIS, "Standard Deviation (sqrt(v))"); range.SetTitleText("Lognormal Parameters", true); range.ShowObs(false); range.SetGraphAreaBackgroundColor( WHITE ); range.SetAspectRatio(1); range.DrawUseDataCoordinates(); range.DrawSetPenColor(BLUE); /* draw horizontal lines */ do i = 1 to ncol(sigma); y0 = repeat(sigma[i], 1, ncol(x)); h = N2LN( x` || y0` ); range.DrawLine(h[,1], h[,2]); end; range.DrawSetPenColor(RED); /* draw vertical lines */ do i = 1 to ncol(mu); x0 = repeat(mu[i], ncol(y)); v = N2LN( x0 || y` ); range.DrawLine(v[,1], v[,2]); end;