Q-values for False Discovey Rate of resels. Usage: qval = SurfStatQ( slm [, mask] ); slm.t = 1 x v vector of test statistics, v=#vertices. slm.df = degrees of freedom. slm.dfs = 1 x v vector of optional effective degrees of freedom. slm.k = #variates. mask = 1 x v logical vector, 1=inside, 0=outside, = ones(1,v), i.e. the whole surface, by default. The following are optional: slm.resl = e x k matrix of sum over observations of squares of differences of normalized residuals along each edge. slm.tri = t x 3 matrix of triangle indices, 1-based, t=#triangles. or slm.lat = 3D logical array, 1=in, 0=out. qval.Q = 1 x v vector of Q-values. qval.mask = copy of mask.
0001 function qval = SurfStatQ( slm, mask ); 0002 0003 %Q-values for False Discovey Rate of resels. 0004 % 0005 % Usage: qval = SurfStatQ( slm [, mask] ); 0006 % 0007 % slm.t = 1 x v vector of test statistics, v=#vertices. 0008 % slm.df = degrees of freedom. 0009 % slm.dfs = 1 x v vector of optional effective degrees of freedom. 0010 % slm.k = #variates. 0011 % mask = 1 x v logical vector, 1=inside, 0=outside, 0012 % = ones(1,v), i.e. the whole surface, by default. 0013 % The following are optional: 0014 % slm.resl = e x k matrix of sum over observations of squares of 0015 % differences of normalized residuals along each edge. 0016 % slm.tri = t x 3 matrix of triangle indices, 1-based, t=#triangles. 0017 % or 0018 % slm.lat = 3D logical array, 1=in, 0=out. 0019 % 0020 % qval.Q = 1 x v vector of Q-values. 0021 % qval.mask = copy of mask. 0022 0023 [l,v]=size(slm.t); 0024 if nargin<2 | isempty(mask) 0025 mask=logical(ones(1,v)); 0026 end 0027 0028 df=zeros(2); 0029 ndf=length(slm.df); 0030 df(1,1:ndf)=slm.df; 0031 df(2,1:2)=slm.df(ndf); 0032 if isfield(slm,'dfs') 0033 df(1,ndf)=mean(slm.dfs(mask>0)); 0034 end 0035 0036 if isfield(slm,'du') 0037 [resels,reselspvert]=SurfStatResels(slm,mask); 0038 else 0039 reselspvert=ones(1,v); 0040 end 0041 reselspvert=reselspvert(mask); 0042 0043 P_val=stat_threshold(0,1,0,df,[10 slm.t(1,mask)],[],[],[],slm.k,[],[],0); 0044 P_val=P_val(2:length(P_val)); 0045 np=length(P_val); 0046 [P_sort, index]=sort(P_val); 0047 r_sort=reselspvert(index); 0048 c_sort=cumsum(r_sort); 0049 P_sort=P_sort./(c_sort+(c_sort<=0)).*(c_sort>0)*sum(r_sort); 0050 m=1; 0051 Q_sort=zeros(1,np); 0052 for i=np:-1:1 0053 if P_sort(i)<m 0054 m=P_sort(i); 0055 end 0056 Q_sort(i)=m; 0057 end 0058 Q=zeros(1,np); 0059 Q(index)=Q_sort; 0060 qval.Q=ones(1,size(mask,2)); 0061 qval.Q(mask)=Q; 0062 qval.mask=mask; 0063 0064 return 0065 end 0066