|
0
|
1 function [COMB_MEAN,COMB_VARIANCE]=get_mean_variance_seg(gene_expression_1,gene_expression_2,region_counts_1,region_counts_2,variance_function_parametric_1, variance_function_parametric_2)
|
|
|
2
|
|
|
3
|
|
|
4 COMB_READS_PER_EXON=[region_counts_1;region_counts_2];
|
|
|
5 GENE_EXPRESSION=[gene_expression_1';gene_expression_2'];
|
|
|
6 IX_SAMPLE1=1:length(gene_expression_1);
|
|
|
7 IX_SAMPLE2=(1+length(gene_expression_1)):(length(gene_expression_1)+length(gene_expression_2));
|
|
|
8
|
|
|
9 COMB_VARIANCE=[];
|
|
|
10 COMB_MEAN=[];
|
|
|
11 for i=1:size(region_counts_1,2)
|
|
|
12 INTEN1=region_counts_1(:,i);
|
|
|
13 INTEN2=region_counts_2(:,i);
|
|
|
14 INTENSITY=[INTEN1;INTEN2];
|
|
|
15
|
|
|
16 SR_VECT=GENE_EXPRESSION;
|
|
|
17
|
|
|
18 %Are there any counts at all in the region?
|
|
|
19 if sum(INTENSITY)>0
|
|
|
20 %Compute the means under the null hypothesis
|
|
|
21 Q=(INTENSITY./SR_VECT)/sum(SR_VECT>0);
|
|
|
22
|
|
|
23 if sum(isnan(SR_VECT))
|
|
|
24 MEAN1=0;
|
|
|
25 MEAN2=0;
|
|
|
26 else
|
|
|
27 MEAN1=mean(sum(Q)*SR_VECT(IX_SAMPLE1));
|
|
|
28 MEAN2=mean(sum(Q)*SR_VECT(IX_SAMPLE2));
|
|
|
29 end
|
|
|
30
|
|
|
31 COMB_MEAN=[COMB_MEAN,[MEAN1;MEAN2]];
|
|
|
32 else
|
|
|
33 COMB_MEAN=[COMB_MEAN,[0;0]];
|
|
|
34 end
|
|
|
35
|
|
|
36 end
|
|
|
37
|
|
|
38 VARIANCE1= predict_variance(COMB_MEAN(1,:)',variance_function_parametric_1)';
|
|
|
39 VARIANCE2= predict_variance(COMB_MEAN(2,:)',variance_function_parametric_2)';
|
|
|
40
|
|
|
41 COMB_VARIANCE=[VARIANCE1;VARIANCE2]; |
