Hello,
I've been trying to write a script that will perform the FDR adjustment at the end of a process where many pvalues are calculated. Several columns of Pvalues are generated and I want to do the FDR adjustment on each columns.
The example I have has 3 columns on which it should calculate the FDR. I recycled a piece of code from the add-in developed by @John_Sall
The script I have has a For Loop to do the job and works well on the first column but the For Loop stops prematurely without processing the other columns and does not throw an error.
I'm at a loss.
Any suggestions would be greatly appreciated.
Thanks
Sebastien
Names Default To Here( 1 );
Clear Log();
New Table( "subset",
Add Rows( 88 ),
New Column( "unique_subject_id_variables",
Character,
"Nominal",
Set Values(
{"ELN-19-0aa3-5-2 # spinal cord", "ELN-19-0aa3-5-3 # spinal cord",
"ELN-19-0aa3-5-4 # spinal cord", "ELN-19-0aa3-5-5 # spinal cord",
"ELN-19-0aa3-5-6 # spinal cord", "ELN-19-0aa3-2-1 # spinal cord",
"ELN-19-0aa3-2-2 # spinal cord", "ELN-19-0aa3-2-3 # spinal cord",
"ELN-19-0aa3-2-4 # spinal cord", "ELN-19-0aa3-2-5 # spinal cord",
"ELN-19-0aa3-2-6 # spinal cord", "ELN-19-0aa3-7-1 # spinal cord",
"ELN-19-0aa3-7-2 # spinal cord", "ELN-19-0aa3-7-3 # spinal cord",
"ELN-19-0aa3-7-4 # spinal cord", "ELN-19-0aa3-7-5 # spinal cord",
"ELN-19-0aa3-7-6 # spinal cord", "ELN-19-0aa3-4-1 # spinal cord",
"ELN-19-0aa3-4-3 # spinal cord", "ELN-19-0aa3-4-5 # spinal cord",
"ELN-19-0aa3-4-6 # spinal cord", "ELN-19-0aa3-3-1 # spinal cord",
"ELN-19-0aa3-3-2 # spinal cord", "ELN-19-0aa3-3-3 # spinal cord",
"ELN-19-0aa3-3-4 # spinal cord", "ELN-19-0aa3-3-5 # spinal cord",
"ELN-19-0aa3-3-6 # spinal cord", "ELN-19-0aa3-6-2 # spinal cord",
"ELN-19-0aa3-6-3 # spinal cord", "ELN-19-0aa3-6-4 # spinal cord",
"ELN-19-0aa3-6-5 # spinal cord", "ELN-19-0aa3-6-6 # spinal cord",
"ELN-19-0aa3-8-1 # spinal cord", "ELN-19-0aa3-8-2 # spinal cord",
"ELN-19-0aa3-8-3 # spinal cord", "ELN-19-0aa3-8-4 # spinal cord",
"ELN-19-0aa3-8-5 # spinal cord", "ELN-19-0aa3-8-6 # spinal cord",
"ELN-19-0aa3-1-1 # spinal cord", "ELN-19-0aa3-1-2 # spinal cord",
"ELN-19-0aa3-1-3 # spinal cord", "ELN-19-0aa3-1-4 # spinal cord",
"ELN-19-0aa3-1-5 # spinal cord", "ELN-19-0aa3-1-6 # spinal cord",
"ELN-19-0aa3-5-2 # brain:cerebellum",
"ELN-19-0aa3-5-3 # brain:cerebellum",
"ELN-19-0aa3-5-4 # brain:cerebellum",
"ELN-19-0aa3-5-5 # brain:cerebellum",
"ELN-19-0aa3-5-6 # brain:cerebellum",
"ELN-19-0aa3-2-1 # brain:cerebellum",
"ELN-19-0aa3-2-2 # brain:cerebellum",
"ELN-19-0aa3-2-3 # brain:cerebellum",
"ELN-19-0aa3-2-4 # brain:cerebellum",
"ELN-19-0aa3-2-5 # brain:cerebellum",
"ELN-19-0aa3-2-6 # brain:cerebellum",
"ELN-19-0aa3-7-1 # brain:cerebellum",
"ELN-19-0aa3-7-2 # brain:cerebellum",
"ELN-19-0aa3-7-3 # brain:cerebellum",
"ELN-19-0aa3-7-4 # brain:cerebellum",
"ELN-19-0aa3-7-5 # brain:cerebellum",
"ELN-19-0aa3-7-6 # brain:cerebellum",
"ELN-19-0aa3-4-1 # brain:cerebellum",
"ELN-19-0aa3-4-3 # brain:cerebellum",
"ELN-19-0aa3-4-5 # brain:cerebellum",
"ELN-19-0aa3-4-6 # brain:cerebellum",
"ELN-19-0aa3-3-1 # brain:cerebellum",
"ELN-19-0aa3-3-2 # brain:cerebellum",
"ELN-19-0aa3-3-3 # brain:cerebellum",
"ELN-19-0aa3-3-4 # brain:cerebellum",
"ELN-19-0aa3-3-5 # brain:cerebellum",
"ELN-19-0aa3-3-6 # brain:cerebellum",
"ELN-19-0aa3-6-2 # brain:cerebellum",
"ELN-19-0aa3-6-3 # brain:cerebellum",
"ELN-19-0aa3-6-4 # brain:cerebellum",
"ELN-19-0aa3-6-5 # brain:cerebellum",
"ELN-19-0aa3-6-6 # brain:cerebellum",
"ELN-19-0aa3-8-1 # brain:cerebellum",
"ELN-19-0aa3-8-2 # brain:cerebellum",
"ELN-19-0aa3-8-3 # brain:cerebellum",
"ELN-19-0aa3-8-4 # brain:cerebellum",
"ELN-19-0aa3-8-5 # brain:cerebellum",
"ELN-19-0aa3-8-6 # brain:cerebellum",
"ELN-19-0aa3-1-1 # brain:cerebellum",
"ELN-19-0aa3-1-2 # brain:cerebellum",
"ELN-19-0aa3-1-3 # brain:cerebellum",
"ELN-19-0aa3-1-4 # brain:cerebellum",
"ELN-19-0aa3-1-5 # brain:cerebellum",
"ELN-19-0aa3-1-6 # brain:cerebellum"}
)
),
New Column( "XYZ % (Norm to Ctrl Gene/Ctrl Grp)_One-sample Prob < t",
Numeric,
"Continuous",
Format( "PValue", 12 ),
Set Values(
[7.00762374713579e-14, 0.0000039742284681085, 0.0120375286563227,
0.0000030636932283006, 0.000252596705471054, 0.958310761173214,
0.932157675214062, 0.606960214034368, 0.99665564502601,
0.0483142647886624, 0.0322337116111752, 0.0000018559647703968,
0.0000000046336081797, 9.85490880976054e-11, 6.96443993397079e-12,
0.0000000019441394764, 7.91203443124632e-13, 0.366497859114321,
0.0000000169366370589, 1.02682759861988e-10, 0.0000006065799989169,
0.103055578812252, 0.995764311763487, 0.980132511993954,
0.995068373759193, 0.984517082520774, 0.0000002181043490997,
0.0000000016031797658, 0.0000000509327294227, 0.0000000087075814375,
0.0000000016497055961, 0.0000000180355378188, 0.0000001520688377862,
5.06777962824243e-11, 3.10032398318472e-10, 0.794729459209574,
0.0000056729575662837, 0.0000033713227982584, 0.920325866551036,
0.579492094647097, 0.947324113072484, 0.00273886144350255,
0.854234625769266, 0.371804857786718, 2.76816294705043e-10,
0.0127493463140254, 0.0787726866182871, 0.00764416297990253,
0.0233627185434858, 0.934031966765819, 0.629482606177694,
0.293765025212661, 0.955254915691655, 0.653020589585128,
0.0361360173941131, 0.409527085611957, 0.753918776098066,
0.992080323172834, 0.906606371236632, 0.987744089041033,
0.0320280732140181, 0.261372779660272, 0.949991060666457,
0.852193284106716, 0.00639756539088631, 0.00311870496140871,
0.944030249515997, 0.918371540034979, 0.23247341672968,
0.445127833217757, 0.123453674637598, 0.164053482001082,
0.798956563445221, 0.0000936772910432205, 0.999345457691391,
0.247703769775099, 0.980401706236753, 0.00246164853294895,
0.0271405807203969, 0.0378473630501983, 0.839608145623234,
0.542075423478042, 0.604346193224194, 0.986203075062795,
0.13778057894564, 0.00761927477682632, 0.964520943546186,
0.230239098364996]
)
),
New Column( "XYZ % (Norm to Ctrl Gene/Ctrl Grp)_One sample Prob > |t|",
Numeric,
"Continuous",
Format( "PValue", 12 ),
Set Values(
[1.40152474942716e-13, 0.0000079484569362171, 0.0240750573126453,
0.0000061273864566012, 0.000505193410942107, 0.0833784776535725,
0.135684649571877, 0.786079571931264, 0.00668870994797958,
0.0966285295773248, 0.0644674232223504, 0.0000037119295407935,
0.0000000092672163595, 1.97098176195211e-10, 1.39288798679416e-11,
0.0000000038882789528, 1.58240688624926e-12, 0.732995718228642,
0.0000000338732741179, 2.05365519723976e-10, 0.0000012131599978339,
0.206111157624505, 0.00847137647302538, 0.0397349760120927,
0.00986325248161459, 0.0309658349584516, 0.0000004362086981994,
0.0000000032063595316, 0.0000001018654588455, 0.000000017415162875,
0.0000000032994111922, 0.0000000360710756376, 0.0000003041376755724,
1.01355592564849e-10, 0.0000000006200647966, 0.410541081580852,
0.0000113459151325673, 0.0000067426455965167, 0.159348266897927,
0.841015810705806, 0.105351773855033, 0.00547772288700511,
0.291530748461467, 0.743609715573435, 0.0000000005536325894,
0.0254986926280507, 0.157545373236574, 0.0152883259598051,
0.0467254370869716, 0.131936066468362, 0.741034787644612,
0.587530050425323, 0.0894901686166905, 0.693958820829744,
0.0722720347882261, 0.819054171223915, 0.492162447803867,
0.0158393536543313, 0.186787257526736, 0.0245118219179341,
0.0640561464280361, 0.522745559320544, 0.100017878667087,
0.295613431786568, 0.0127951307817726, 0.00623740992281741,
0.111939500968007, 0.163256919930042, 0.464946833459361,
0.890255666435514, 0.246907349275197, 0.328106964002164,
0.402086873109557, 0.000187354582086441, 0.00130908461721828,
0.495407539550197, 0.0391965875264939, 0.00492329706589789,
0.0542811614407937, 0.0756947261003966, 0.320783708753531,
0.915849153043917, 0.791307613551612, 0.0275938498744102,
0.27556115789128, 0.0152385495536526, 0.0709581129076288,
0.460478196729993]
)
),
New Column( "XYZ % (Norm to Ctrl Gene/Ctrl Grp)_One sample Prob > t",
Numeric,
"Continuous",
Format( "PValue", 12 ),
Set Values(
[0.99999999999993, 0.999996025771532, 0.987962471343678,
0.999996936306772, 0.999747403294529, 0.0416892388267862,
0.0678423247859383, 0.393039785965632, 0.00334435497398979,
0.951685735211338, 0.967766288388825, 0.99999814403523,
0.999999995366392, 0.999999999901451, 0.999999999993036,
0.999999998055861, 0.999999999999209, 0.633502140885679,
0.999999983063363, 0.999999999897317, 0.999999393420001,
0.896944421187748, 0.00423568823651269, 0.0198674880060463,
0.0049316262408073, 0.0154829174792258, 0.999999781895651,
0.99999999839682, 0.999999949067271, 0.999999991292419,
0.999999998350294, 0.999999981964462, 0.999999847931162,
0.999999999949322, 0.999999999689968, 0.205270540790426,
0.999994327042434, 0.999996628677202, 0.0796741334489637,
0.420507905352903, 0.0526758869275166, 0.997261138556498,
0.145765374230734, 0.628195142213282, 0.999999999723184,
0.987250653685975, 0.921227313381713, 0.992355837020097,
0.976637281456514, 0.0659680332341809, 0.370517393822306,
0.706234974787339, 0.0447450843083453, 0.346979410414872,
0.963863982605887, 0.590472914388043, 0.246081223901934,
0.00791967682716566, 0.0933936287633682, 0.0122559109589671,
0.967971926785982, 0.738627220339728, 0.0500089393335434,
0.147806715893284, 0.993602434609114, 0.996881295038591,
0.0559697504840033, 0.0816284599650208, 0.76752658327032,
0.554872166782243, 0.876546325362402, 0.835946517998918,
0.201043436554779, 0.999906322708957, 0.000654542308609138,
0.752296230224901, 0.019598293763247, 0.997538351467051,
0.972859419279603, 0.962152636949802, 0.160391854376766,
0.457924576521958, 0.395653806775806, 0.0137969249372051,
0.86221942105436, 0.992380725223174, 0.0354790564538144,
0.769760901635004]
)
)
);
dx=data table("subset");
c = "Prob ";
col = dx << get column names( string );
show (col);
Probcols = {};
nc = N Items( col );
For( i = 1, i <= nc, i++,
If( Contains( col[i], c ),show(col[i]);
insert into(Probcols,col[i])
));
show(Probcols); //list of Prob columns
pvalnbr=N items(Probcols);
show(pvalnbr);
for(i=1,i<=pvalnbr,i++,
show(Probcols[i]);//there is 3 columns to perform FDR calculation on but the loop stops after 1 round
pvalues = Column(Probcols[i])<<Get Values;
if (min(pvalues)<0 | max(pvalues)>1,dialog("PValues must be between 0 and 1");throw());
// ranks the p-values and applies Benjamini-Hochberg adjustment to them - Code from FDR add-in developped by John Sall@JMP
index = rankIndex(pvalues); // removes missing values
n = nrow(index);
ordPValue = pvalues[index];
adjPValue = J(n,1,.);
for(i=1,i<=n,i++,
j = n-i+1;
ratio = n/(n-i+1);
adjPValue[j] = if (i==1,ordPValue[j],min(ratio*ordPValue[j],adjPValue[j+1]));
);
rowPValues = J(NRow(dx),1,.);
rowPValues[index] = adjPValue;
newcol=dx<<New Column("FDR Adj PValue",Numeric,Values(rowPValues));
//newcol<<set name(Probcols[i]||" FDR Adj PValue") // need to also figure out why it does not rename - likely related to the reason why it coming out of the loop prematurely
);