Re: Computing Wilcoxon signed rank test p-value
- From: "Gene" <Nospam54@xxxxxxx>
- Date: 18 Jan 2006 07:36:49 -0800
Hollander & Wolfe has a good description. The only tricky part is
finding the pattern of ties in the data. Here's my little Matlab
implementation. The ties.m program is my own and it's slow, but it
works .
function [pvalue,W]=Wilcoxsignrank(X,Y)
% Wilcoxon signed rank test
% [pvalue,W]=Wilcoxsignrank(X,Y)
% Tests the null hypothesis that X & Y have the same pdf.
% Input: X,Y two paired samples
% Output: pvalue:2-sided p value for large sample approximation N(0,1)
distribution
% W=Wilcoxon signed rank statistic, sum of the positive ranks.
% Calls Gallagher's ties.m (included below) to find the pattern of
ties.
% Reference: Hollander & Wolfe (1999) Nonparametric statistical methods
% Larsen & Marx (2001), page 692-695
% Written by E. Gallagher;
X=X(:);Y=Y(:);
n=length(X);m=length(Y);
if n~=m
error('X & Y must be the same length')
end
% Rank the abs(X-Y) values from smallest to largest, assigning average
ranks to ties.
Diff=X-Y;Psi=Diff>0;
[T,R,ind]=ties(abs(Diff));T=T';R=R'; % calls Gallagher's ties.m
W=sum(Psi.*R); % Hollander & Wolfe, Equation 3.3.
% Large sample approximation;% Hollander & Wolfe p. 108
EoW=(n*(n+1))/4;
% Calculate the variance of W, without ties and with ties.
if isempty(T) % Size of tied groups from ties.m
VaroW=(n*(n+1)*(2*n+1))/24;
else
VaroW=( n*(n+1)*(2*n+1) -0.5*(sum((T-1).*T.*(T+1))) ) / 24;
end
Tstar=(W-EoW)/sqrt(VaroW); % Without ties, tends to an asymptotic
N(0,1) distribution.
% Find the 2-tailedprobability of Wstar from the standard normal
distributioin
pvalue=erfc(abs(Tstar)/sqrt(2));
% Note that the exact p values are tabulated, and an exact test, even
in the presence of ties
% can be performed, see pp. 113-116 in Hollander & Wolfe.
function [T,R,ind]=ties(A)
% format: [T,R,ind]=ties(A)
% a function to return a row vector of tied groups, T,
% Ranks R (including average ranks) and indices of tied elements
% needed to calculate variance of S using Kendall's
% variance formula & Spearman's r.
% input: A is a row or column vector
% T: a row vector containing number of members of tied groups
% T=0 if there are no tied groups
% sum(T) is equal to the number of tied elements.
% each element of T equals the number in each tied group
% tied groups are sorted in ascending order.
% Examples: A=[1 2 3];[T,R,ind]=ties(A)=> T=0,R=[1 2 3],ind=[]
% A=[1 2 3 1]; T=2,R=[1.5 3 4 1.5],ind=[1 4]
% A=[2 1 2 3 1 2]; T=[2 3],R=[4 1.5 4 6 1.5 4],
% ind=[5 2 3 1 6]
% A=[2 1 2 3 3 1 2]; T=[2 3 2],R=[4 1.5 4 6.5 6.5 1.5 4]
% ind=[6 2 3 1 7 4 5]
% R (Row vec)=numerical rankings of A with ave. ranks for ties
% ind: indices of tied elements, sorted by rank; sorted tied
elements=A(ind);
% ties.m is used in Kendall.m as T=ties(A), and Spear.m
% written by E. Gallagher, Environmental Sciences Program
% UMASS/Boston, Email: Eugene.Gallagher@xxxxxxx
% written: 6/16/93, revised 6/17/93
[r,c]=size(A);
if r>c
A=A'; % change to row vector
end
[Asort,k]=sort(A);
iota=1:length(A);iota=iota';
R(k)=iota;
index=[k' iota];
ind=[];
CDA=[~diff(Asort) 0];
min1=min(find(CDA==1));
if isempty(min1)
T=0;
return
end
i=0;
[rw,cl]=size(CDA);
T=zeros(size(rw,cl));
while ~isempty(min1)
min0=min(find(CDA==0));
if min0<min1
CDA(min0:min1-1)=[];
index(min0:min1-1,:)=[];
else
i=i+1;
T(i)=min0-min1+1;
CDA(min1:min0)=[];
ind=[ind index(min1:min0,1)'];
R(1,index(min1:min0))=ones(1,T(i))*sum(index(min1:min0,2))/T(i);
index(min1:min0,:)=[];
end
min1=min(find(CDA==1));
end
T(find(T==0))=[];
.
- Follow-Ups:
- Re: Computing Wilcoxon signed rank test p-value
- From: Pawel Kusmierek
- Re: Computing Wilcoxon signed rank test p-value
- References:
- Computing Wilcoxon signed rank test p-value
- From: Pawel Kusmierek
- Computing Wilcoxon signed rank test p-value
- Prev by Date: Re: experts in small probabilities
- Next by Date: Re: a tricky question on probabilities
- Previous by thread: Computing Wilcoxon signed rank test p-value
- Next by thread: Re: Computing Wilcoxon signed rank test p-value
- Index(es):
Relevant Pages
|
