MATLAB Answers

Dr. Seis
0

Plotting lognormal data on this type of graph?

Asked by Dr. Seis
on 30 Aug 2012
Latest activity Answered by Kirsten Koehler on 15 May 2019
I would like to plot data on a graph similar to that illustrated in the figure below (except it is upside-down from my convention - what they call "99" on the Cumulative probability axis I would call "1" and vice versa for "1" and "99", respectively):
Notice that the divisions from 99 to 50 is a mirror image of the divisions from 50 to 1.
Any ideas? Is there something built in Matlab that does this type of thing already? I am trying to determine P10, P50, P90 information from data distributed lognormally.

  0 Comments

Sign in to comment.

2 Answers

Answer by Dr. Seis
on 31 Aug 2012
Edited by Dr. Seis
on 2 Sep 2012
 Accepted Answer

The key ended up being norminv. Here is an example output of the program I made (see below for probit function):
>> probit([10,90],[500,2000]);
function probit(PP,xx)
if min(PP) >= 1
PP = PP/100;
end
% Define x- and y-axis ranges
x = [.01,100000];
y = [1,1.3,2,5,10,20,23,30,40,50,60,70,77,80,90,95,98,98.7,99];
% Define y-axis in terms of norminv
y_prime = norminv(y/100,0,1);
figure;
% Create probit axis
for i = 1 : numel(y)
if any(y(i) == [1,10,50,90,99])
semilogx(x,[1,1]*y_prime(i),'k-','LineWidth',2);
elseif any(y(i) == [1.3,23,77,98.7])
semilogx(x,[1,1]*y_prime(i),'k--');
else
semilogx(x,[1,1]*y_prime(i),'k:');
end
if (i == 1); hold on; end
end
% Invert input data for slope and y-intercept
GG = ones(numel(PP),2);
for i = 1 : numel(PP)
GG(i,1) = log10(xx(i));
end
temp = GG\reshape(norminv(PP,0,1),numel(PP),1);
slope = temp(1);
yint = temp(2);
% Define best fit line in Least Squares Sense
temp = (y_prime - yint)/slope;
newxx = 10.^temp;
% Plot best fit lines and original/defined points
semilogx(newxx,y_prime,'b-','LineWidth',2);
h = zeros(1,numel(y)); cellstr1={''};
temp = numel(y):-1:1;
for i = 1:numel(y)
h(i) = semilogx(newxx(temp(i)),y_prime(temp(i)),'ro','MarkerSize',15);
cellstr1(i) = cellstr(sprintf('P%4.1f = %10.4f',y(temp(i)),newxx(temp(i))));
end
temp = norminv(PP,0,1);
for i = 1 : numel(PP)
h(numel(y)+1)=semilogx(xx(i),temp(i),'gs','MarkerSize',15);
end
cellstr1(numel(y)+1) = {'Original Points'};
temp = floor(log10(min(newxx))):ceil(log10(max(newxx)));
for i = 1 : numel(temp)
semilogx(10^temp(i)*[1,1], y_prime([1,end]),'k-','LineWidth',2);
end
hold off;
xlim([10^floor(log10(min(newxx))),10^ceil(log10(max(newxx)))]);
ylim([norminv(.01,0,1),norminv(.9901,0,1)]);
set(gca,'YTick',y_prime);
set(gca,'YTickLabel',num2cell(y));
set(gca,'XGrid','on');
legend(h,cellstr1,'Location','NorthWest'); set(gca,'FontSize',14);

  0 Comments

Sign in to comment.


Answer by Kirsten Koehler on 15 May 2019

Thank you for this! I've been taught this is called a log-probability plot, so i'm commenting to add these words because it took me hours to find your code. In my field it's common to swap the axes, but a simple
view([90 -90])
takes care of that.

  0 Comments

Sign in to comment.