analyzeSampling Compares the significance of change in flux between two conditions with the significance of change in gene expression Tex a vector of t-scores for the change in gene expression for each reaction. This score could be the Student t between the two conditions, or you can calculate it from a p-value (by computing the inverse of the so called error function). If you choose the second alternative you should be aware that the transcripts that increased in expression level should have positive values and those who decreased in expression level should have negative values (the p-values only tell you if the fluxes changed or not but not in which direction) df the degrees of freedom in the t-test solutionsA random solutions for the reference condition (as generated by randomSampling) solutionsB random solutions for the test condition (as generated by randomSampling) printResults prints the most significant reactions in each category (opt, default false) scores a Nx3 column matrix with the probabilities of a reaction: 1) changing both in flux and expression in the same direction 2) changing in expression but not in flux 3) changing in flux but not in expression or changing in opposed directions in flux and expression. Usage: scores=analyzeSampling(Tex, df, solutionsA, solutionsB, printResults) Rasmus Agren, 2012-03-19
0001 function scores=analyzeSampling(Tex, df, solutionsA, solutionsB, printResults) 0002 % analyzeSampling 0003 % Compares the significance of change in flux between two conditions with 0004 % the significance of change in gene expression 0005 % 0006 % Tex a vector of t-scores for the change in gene expression 0007 % for each reaction. This score could be the Student t 0008 % between the two conditions, or you can calculate it from 0009 % a p-value (by computing the inverse of the so called error 0010 % function). If you choose the second alternative you should 0011 % be aware that the transcripts that increased in expression 0012 % level should have positive values and those who decreased 0013 % in expression level should have negative values (the 0014 % p-values only tell you if the fluxes changed or not but 0015 % not in which direction) 0016 % df the degrees of freedom in the t-test 0017 % solutionsA random solutions for the reference condition (as 0018 % generated by randomSampling) 0019 % solutionsB random solutions for the test condition (as generated 0020 % by randomSampling) 0021 % printResults prints the most significant reactions in each category 0022 % (opt, default false) 0023 % 0024 % scores a Nx3 column matrix with the probabilities of a reaction: 0025 % 1) changing both in flux and expression in the same direction 0026 % 2) changing in expression but not in flux 0027 % 3) changing in flux but not in expression or changing 0028 % in opposed directions in flux and expression. 0029 % 0030 % Usage: scores=analyzeSampling(Tex, df, solutionsA, solutionsB, printResults) 0031 % 0032 % Rasmus Agren, 2012-03-19 0033 % 0034 0035 if nargin<5 0036 printResults=false; 0037 end 0038 0039 nRxns=numel(Tex); 0040 pM=zeros(nRxns,1); 0041 pH=zeros(nRxns,1); 0042 pR=zeros(nRxns,1); 0043 0044 %Check that the number of reactions is the same in both expression and flux 0045 if nRxns~=size(solutionsA,1) 0046 throw(MException('','The number of reactions must be the same in Tex as in solutionsA')); 0047 end 0048 0049 %Get the Z-score and mean for the solutions 0050 mA=mean(solutionsA,2); 0051 mB=mean(solutionsB,2); 0052 Zf=getFluxZ(solutionsA, solutionsB); 0053 0054 %Clear up the tex if there are elements that are NaN or +/- Inf. 0055 I=isnan(Tex) | isinf(Tex); 0056 if any(I) 0057 fprintf('WARNING: There are t-scores that are NaN or +/- Inf. These values are changed to 0.0\n'); 0058 end 0059 Tex(I)=0; 0060 0061 for i=1:nRxns 0062 %Check how the flux has changed. The means are needed because to 0063 %differentiate between a positive flux changing to a smaller flux and a 0064 %negative flux changing to a more negative flux (which is a larger 0065 %flux) 0066 I=mB(i)/mA(i); 0067 if I<0 0068 pM(i)=erf(abs(Zf(i))); 0069 pH(i)=(1-pM(i))*(2*tcdf(abs(Tex(i)),df)-1); 0070 pR(i)=0; 0071 else 0072 if mB(i)<0 0073 Zf(i)=Zf(i)*-1; 0074 end 0075 end 0076 0077 I=Zf(i)/Tex(i); 0078 if I<0 0079 pM(i)=erf(abs(Zf(i))); 0080 pH(i)=(1-pM(i))*(2*tcdf(abs(Tex(i)),df)-1); 0081 pR(i)=0; 0082 else 0083 pR(i)=erf(abs(Zf(i)))*(2*tcdf(abs(Tex(i)),df)-1); 0084 pH(i)=(2*tcdf(abs(Tex(i)),df)-1)*(1-erf(abs(Zf(i)))); 0085 pM(i)=erf(abs(Zf(i)))*(1-(2*tcdf(abs(Tex(i)),df)-1)); 0086 end 0087 end 0088 0089 scores=[pR pH pM]; 0090 0091 if printResults==true 0092 fprintf('TOP SCORING REACTIONS\n\n'); 0093 %The top 10 hits in the first category 0094 [I J]=sort(pR,'descend'); 0095 fprintf('Reactions which change both in flux and expression in the same direction\nReaction\tProbability\n'); 0096 for i=1:10 0097 fprintf([num2str(J(i)) '\t' num2str(I(i)) '\n']); 0098 end 0099 0100 %The top 10 hits in the first category 0101 [I J]=sort(pH,'descend'); 0102 fprintf('\nReactions which change in expression but not in flux\nReaction\tProbability\n'); 0103 for i=1:10 0104 fprintf([num2str(J(i)) '\t' num2str(I(i)) '\n']); 0105 end 0106 0107 %The top 10 hits in the first category 0108 [I J]=sort(pM,'descend'); 0109 fprintf('\nReactions which change in flux but not in expression, or in opposed directions in flux and expression\nReaction\tProbability\n'); 0110 for i=1:10 0111 fprintf([num2str(J(i)) '\t' num2str(I(i)) '\n']); 0112 end 0113 end 0114 end