(*:Mathematica:: Version 7*)
(*:Context: "swirl`"         *)
(*:Swirl-and-Ripple       *)
(*:Version                1*)
(*Title: Supporting Functions*)
(*:Author: Stuart G. Baker,  2010                   *)



BeginPackage["swirlsup`"]

RanSample::usage="RanSample[data,split,showprogress]"

StatAUCList::usage="StatAUCList[x0,x1], where x1 and x1 are matrices"
StatAUC::usage="StatAUC[xvec0,xvec1]"

SelectDistanceMeasure::usage="SelectDistanceMeasure[{auc1,auc2},delta]"

Var::usage="Var[list]"
VarX::usage="VarX[list] is variance with adjustment for a list of length one"
CDFNormal01::usage="CDFNormal01[x]"
ReplaceZero::usage="ReplaceZero"
SelectRows::usage="SelectRows"
SelectColumns::usage="SelectColumns"



Clear[RanSample, SelectRows, SelectTopGenes,StatAUClist,StatAUC, 
     SelectDistanceMeasure,VarX,Var,CDFNormal01,ReplaceZero]
	  

Begin["Private`"]
														
     



(*-------------Select distance measure------------------------------------------*)



SelectDistanceMeasure[{auc1_,auc2_},delta_]:=
Module[{len},
If[auc2-auc1 < delta, len=1];
If[auc2-auc1 >=delta, len=2];
Return[len]]



(*------------Random sample for split-------------------------------*)
 
RanSample[{mat0_,mat1_},split_,shownumsplit_]:=
Module[{train0,test0,train1,test1,res,trainsubjects0,trainsubjects1,mat0x,mat1x},
  {train0,test0,trainsubjects0}=ransplitroc[mat0,split];
  {train1,test1,trainsubjects1}=ransplitroc[mat1,split];
  res={{train0,train1},{test0,test1},{trainsubjects0,trainsubjects1}};
Return[res]]



ransplitroc[matx_,split_]:=
 Module[{mat,i,len,range,sample,pairs,
     lenTRAIN,lenTEST,
     pairsTRAIN,pairsTEST,
	  sampleTRAIN,sampleTEST,
	  matTRAIN,matTEST,res},	  
  (*create matrix of subjects by genes from matx which genes by subject*)
  	    mat=Transpose[matx];	
	(*lengths*)
             len=Length[mat];
	    lenTRAIN=Max[1,Round[len split]];
	    range=Range[len];
	 (*random positions in sample*)
	    sample=Table[Random[],{i,1,len}];
	 (*reorder range by random sample *)
	    pairs=Sort@Transpose[{sample,range}];
	(*for training sample TAKE lenTRAIN*)
	     pairsTRAIN=Take[pairs,lenTRAIN];
	 (*for test sample DROP lenTRAIN*) 
	  	  pairsTEST=Drop[pairs,lenTRAIN];
	 (*idenity new ranges*)
	 	 {sampleTRAIN,rangeTRAIN}=Transpose[pairsTRAIN];
	     {sampleTEST,rangeTEST}=Transpose[pairsTEST];
	 (*take values according to new ranges*)
 	 	  matTRAIN=mat[[rangeTRAIN]];
	     matTEST=mat[[rangeTEST]];
	res={Transpose[matTRAIN],Transpose[matTEST],rangeTRAIN};  (*genes by subject*)
  Return[res]]



 (*------------AUC of list based on binormal distribution------------------------------------*)
 
 
 StatAUCList[y0_,y1_]:=MapThread[StatAUC[#1,#2]&, {y0,y1}]


StatAUC[y0_,y1_]:=
 Module[{m0,m1,v0,v1,n0,n1,auc},
  m0=Mean[y0];
  m1=Mean[y1];
  n0=Length[y0];
  n1=Length[y1];
  If[n0 >1 && n1>1,
   v0=ReplaceZero[Var[y0],.01];
   v1=ReplaceZero[Var[y1],.01];
   t=Abs[(m0-m1)/Sqrt[v0+v1]]//N,
   v0=.01;
   v1=.01];
   auc=CDFNormal01[t]; 
  Return[auc]]	


CDFNormal[mu_,sigma_,x_]:=N[(Erf[(x-mu)/(Sqrt[2] sigma)] + 1) / 2]
CDFNormal01[x_]:=CDFNormal[0,1,x]

(*---------------Other Functions----------------------------------*)

(*variance adjusted if only one*)

VarX[y_]:=
Module[{v},
  If[Length[y]>1, v=Var[y], Print["warning length of data=1"]; v=.01];
Return[v]]


 Var[y_]:= Apply[Plus,(y-Mean[y])^2]/(Length[y]-1)//N;


ReplaceZero[x_,new_:1]:=
   x + ZeroSign[x] new   

ZeroSign[x_]:=1-Abs[Sign[x]]


SelectRows[matrix_List?MatrixQ,indices_List?VectorQ] := matrix[[#]]& /@ indices

SelectColumns[matrix_List?MatrixQ, indices_List?VectorQ] :=           
   Transpose @ SelectRows[Transpose[matrix],indices]  



End[] 
EndPackage[]

