全部版块 我的主页
论坛 数据科学与人工智能 数据分析与数据科学 MATLAB等数学软件专版
5045 3
2015-06-03
悬赏 15 个论坛币 未解决
现在有从论坛上下载的moran 的代码,如下:求问这个代码如何修正,运行总是有问题。这个代码只能做单变量的moran 吗?

function[I,I_standard,Z_I,Z_alpha,result]=moran(X,W,alpha)

%%X为列向量,W为权重矩阵

%alpha为显著性水平

n = length(X);

Xmean = mean(X);

XminusX_mean = X - Xmean;

S0 = sum(sum(W));


%一、计算,moran'I

%1W未标准化,为0-1矩阵

I =XminusX_mean'*W*XminusX_mean*n/(XminusX_mean'*XminusX_mean*S0);

%2W标准化,为(0,1)矩阵

W_rowSum = sum(W');

W_rowSum_temp = W_rowSum'*ones(1,n);

W_standard = W./W_rowSum_temp;

I_standard =XminusX_mean'*W_standard*XminusX_mean/(XminusX_mean'*XminusX_mean);


%二、moran'I值的显著性检验

Z_I = moran_test(I_standard,W_standard);

Z_alpha = norminv(1-alpha,0,1);

%三、画moran散点图

%y = ax + b;线性拟合

%1X标准化,W标准化

   Xstd = std(X,1);%X的标准差

   X_std = (X-Xmean)/Xstd;%X标准化

   W_standard_X_standard = W_standard*X_std;

    %result = myls(X_std,WX_standard);

   result = moranScatterPlot(X_std, W_standard_X_standard );

   title('X标准化,W标准化');


% %2X标准化,W未标准化

%    Xstd = std(X,1);%X的标准差

%    X_std = (X-Xmean)/Xstd;%X标准化,X_std代表横坐标

%    W_X_standard = W*X_std;%WX_standard代表纵坐标

%    % result = myls(X_std,WX_standard);

%   result =  moranScatterPlot(X_std,W_X_standard);

%   title('X标准化,W未标准化');


%3X为未标准化,W标准化

%    W_standard_X = W_standard*X;%WX_standard代表纵坐标

%    % result = myls(X_std,WX_standard);

%    result = moranScatterPlot(X,W_standard_X)

%      title('X为未标准化,W标准化');


% % 4X未标准化,W未标准化

%    WX = W*X;%WX_standard代表纵坐标

%    result = moranScatterPlot(X,WX);

%    title('X未标准化,W未标准化');


Moran I 检验

function [Z_I,Z_alpha] =moran_test(I,W,alpha)

%%I为计算出来的moran'I

%W为计算相应moran'I值的权重矩阵

n = size(W,1);

E_I = -1/(n-1);

S0 = sum(sum(W));

S1 = sum( sum( (W+W').^2 ) )/2;

S2 =sum( sum((W+W')').^2 );

Var_I =(n*n*S1-n*S2+3*S0*S0)/((n*n-1)*S0*S0)-E_I*E_I;

Z_I = (I-E_I)/sqrt(Var_I);

Z_alpha = norminv(1-alpha,0,1);



Moran I plot

function result = moranScatterPlot(X,WX)

   %%X为列向量,W为权重矩阵

   result = myls(X,WX);

   X_lift = min(X) - 0.5*abs(min(X));

   X_right = max(X) + 0.5*max(X);

   xx = linspace(X_lift,X_right,1000);

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

   plot(X,WX,'ro');

   hold on

   plot(xx,yy)

   grid on

   axis equal

end


二维码

扫码加我 拉你入群

请注明:姓名-公司-职位

以便审核进群资格,未注明则拒绝

全部回复
2015-10-17 22:44:23
这个代码确实只能做单变量的哈,我也是这几天在研究,代码有一些小问题,个人觉得
二维码

扫码加我 拉你入群

请注明:姓名-公司-职位

以便审核进群资格,未注明则拒绝

2018-5-20 00:02:23
请问你选择解决了吗,同问题
二维码

扫码加我 拉你入群

请注明:姓名-公司-职位

以便审核进群资格,未注明则拒绝

2018-8-9 18:23:44
同问,你的问题解决了没
二维码

扫码加我 拉你入群

请注明:姓名-公司-职位

以便审核进群资格,未注明则拒绝

相关推荐
栏目导航
热门文章
推荐文章

说点什么

分享

扫码加好友,拉您进群
各岗位、行业、专业交流群