北太天元上用脚本简单实现qqplot的基本功能

北太天元目前还没有实现qqplot画图函数，在很多软件中，这个函数往往属于统计工具箱的一部分。 我先简单实现了一下这个函数， 供大家参考。

qqplot(z)  会画出多个图来，每个图对应z的一列， 每个图的横坐标是标准正态分布的分位数，纵坐标是 z 的每一列数据 分位数, 虚线是线性拟合得到结果。

qqplot(x,y) 假设 x 和 y 是长度一致的向量， 横坐标和纵坐标分别对应 x 和 y 的分位数。

下面是是qqplot 的qqplotEx1.m 的北太天元脚本，含有两个子函数，一个用来实现qqplot, 另一个是

求数据的线性拟合曲线。

在北太天元下运行，得到qqplot图如下

clf

close all

clear all

M=100;N=1;

x=normrnd(0,1,M,N);

y=rand(M,N);

z=[x,y];

figure(1)

qqplot(z);

xlabel('标准正态样本的Quantile');

ylabel('输入样本的Quantile');

title('Quantile-Quantile图');

legend('正态分布数据','均匀分布数据');

%生成两个正态分布样本

x=normrnd(0,1,100,1);

y=normrnd(0.5,2,100,1);

figure(2);

qqplot(x,y);

xlabel('输入样本x的Quantile');

ylabel('输入样本y的Quantile');

title('Quantile-Quantile图');

function qqplot(x,y)

if (nargin == 1)

z = x;

[m,n] = size(z);

t = (double(1:m) -0.5)/m;  %分位数

x = norminv(t);  % norminv 是 标准正态分布的累积分布函数的逆

z = sort(z)

for k=1:n

    plot(x, z(:,k) )

    hold on

    yy = linear_polyfit(x, z(:,k), x);

    plot(x, yy, 'r--')

end

hold off

elseif (nargin == 2)

x = sort(x(:) );

y = sort(y(:) );

length_x = length(x); length_y =  length(y);

if ( length_x ~= length_y )

        error("qqplot(x,y) 的 x 和y 的长度不一样, 需要改成一样");

end

yy = linear_polyfit(x, y, x);

plot(x, y) ;

hold on

plot(x,yy, 'r--')

hold off

else

error("qqplot 目前仅仅处理一个输入参数或者两个输入参数的情况");

end

end

   

function yy =  linear_polyfit(x,y, xx)

p =polyfit(x(:),y(:),1);

yy = p(1)*xx + p(2) ;

end