/usr/share/octave/packages/statistics-1.3.0/normplot.m is in octave-statistics 1.3.0-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 | ## Author: Paul Kienzle <pkienzle@users.sf.net>
## This program is granted to the public domain.
## -*- texinfo -*-
## @deftypefn {Function File} normplot (@var{X})
## Produce normal probability plot for each column of @var{X}.
##
## The line joing the 1st and 3rd quantile is drawn on the
## graph. If the underlying distribution is normal, the
## points will cluster around this line.
##
## Note that this function sets the title, xlabel, ylabel,
## axis, grid, tics and hold properties of the graph. These
## need to be cleared before subsequent graphs using 'clf'.
## @end deftypefn
function normplot(X)
if nargin!=1, print_usage; end
if (rows(X) == 1), X=X(:); end
# Transform data
n = rows(X);
if n<2, error("normplot requires a vector"); end
q = norminv([1:n]'/(n+1));
Y = sort(X);
# Find the line joining the first to the third quartile for each column
q1 = ceil(n/4);
q3 = n-q1+1;
m = (q(q3)-q(q1))./(Y(q3,:)-Y(q1,:));
p = [ m; q(q1)-m.*Y(q1,:) ];
# Plot the lines one at a time. Plot the lines before overlaying the
# normals so that the default label is 'line n'.
if columns(Y)==1,
leg = "+;;";
else
leg = "%d+;Column %d;";
endif
for i=1:columns(Y)
plot(Y(:,i),q,sprintf(leg,i,i)); hold on;
# estimate the mean and standard deviation by linear regression
# [v,dv] = wpolyfit(q,Y(:,i),1)
end
# Overlay the estimated normal lines.
for i=1:columns(Y)
# Use the end points and one point guaranteed to be in the view since
# gnuplot skips any lines whose points are all outside the view.
pts = [Y(1,i);Y(q1,i);Y(end,i)];
plot(pts, polyval(p(:,i),pts), [num2str(i),";;"]);
end
hold off;
# plot labels
title "Normal Probability Plot"
ylabel "% Probability"
xlabel "Data"
# plot grid
t = [0.00001;0.0001;0.001;0.01;0.1;0.3;1;2;5;10;25;50;
75;90;95;98;99;99.7;99.9;99.99;99.999;99.9999;99.99999];
set(gca, "ytick", norminv(t/100), "yticklabel", num2str(t));
grid on
# Set view range with a bit of space around data
miny = min(Y(:)); minq = min(q(1),norminv(0.05));
maxy = max(Y(:)); maxq = max(q(end),norminv(0.95));
yspace = (maxy-miny)*0.05; qspace = (q(end)-q(1))*0.05;
axis ([miny-yspace, maxy+yspace, minq-qspace, maxq+qspace]);
end
|