Dinh/Dinh 2012/NOTES/2012-2-21

From ZhangLabWiki
Jump to navigation Jump to search

48 African Samples - Cleaning[edit]

  • Cleaning with Tech_Rep_9 data: Using requirement of (observed STDEV > 0.1 and p-value < 0.05) I obtained a list of 107,218 bad_cpgs (out of 412,476). The positions are 1-based and on Hg19.
  • Next, I filtered out the sites using this script: File:Getsites.txt (This should be used temporarily, we need a more standardized script). The current UPENN methylMatrix were generated using Dr. Zhang's script so it has 0-based positions.
  • wc -l:
  56092 cleanUPENN46_hg19_methylMatrix.txt
  77573 UPenn46_hg19_methylMatrix_filtered.txt
  • Also use good_cpg_0.1 list (where only the CpGs with p-value > 0.1) to obtain 15,168 CpGs in UPenn data which are high confidence.

Comparing before and after cleaning[edit]

  • Before cleaning:
> A = read.table("UPenn46_hg19_methylMatrix_filtered.txt", header=TRUE, row.names=1)
> B = na.omit(A)
> str(B)
'data.frame':   53548 obs. of  46 variables:
...
> B$STDEV=NULL
> B$min_RD=NULL
> B$mean_RD=NULL
> pca_sites = prcomp(B, scale=TRUE)
> summary(pca_sites)
...                    PC1   PC2    PC3    PC4    PC5     PC6     PC7
Proportion of Variance 0.723 0.0376 0.0125 0.0105 0.00979 0.00869 0.0076
...
> write.table(file="pca_sites.orig.rotations", pca_sites$rotation)
> pca_indiv = prcomp(t(B), scale=TRUE)
> summary(pca_indiv)
...                      PC1    PC2     PC3     PC4     PC5     PC6     PC7
Proportion of Variance   0.205  0.0382  0.0376  0.0336  0.0327  0.0230  0.0216
...
> write.table(file="pca_indiv.orig.x", pca_indiv$x)
  • After cleaning:
> A = read.table("cleanUPENN46_hg19_methylMatrix.txt", header=TRUE, row.names=1)
> C = na.omit(A)
> C$STDEV=NULL
> C$min_RD=NULL
> C$mean_RD=NULL
> str(C)
...
'data.frame':   37135 obs. of  46 variables:
...
> pca_sites = prcomp(C, scale=TRUE)
> summary(pca_sites)
...                    PC1   PC2    PC3    PC4    PC5    PC6     PC7 
Proportion of Variance 0.721 0.0394 0.0125 0.0106 0.0092 0.00861 0.00761
> write.table(file="pca_sites.clean.rotations", pca_sites$rotation)
> pca_indiv = prcomp(t(C), scale=TRUE)
> summary(pca_indiv)
...                     PC1    PC2     PC3     PC4     PC5    PC6     PC7
Proportion of Variance  0.212  0.0378  0.0366  0.0332  0.030  0.0231  0.0217
> write.table(file="pca_indiv.clean.x", pca_indiv$x)

Results:[edit]

 PCA_sites: calculate the correlation of the rotation (weight) of each samples onto CpG sites with QC factor
                   PC1        PC2         PC3
 Original data:    
 %variance         72.3       3.76        1.25
 Pearson's cor(R)  0.77       0.73        0.10
 Clean data:
 %variance         72.1       3.94        1.25
 Pearson's cor(R)  0.77       0.72        0.07
 Clean data (autosomes):
 %variance         71.5       3.82        1.27
 Pearson's cor(R)  0.77       0.73        0.07
 BEST data (autosomes):
 %variance         70.4       4.33        1.37
 Pearson's cor(R)  0.71       0.72        0.10
 PCA_indivs: calculate the correlation of the rotated value of each samples from linear combination of CpG sites with QC factor
                   PC1        PC2         PC3
 Original data:    
 %variance         20.5       3.82        3.76
 Pearson's cor(R)  -0.74      -0.23       -0.45
 Clean data:
 %variance         21.2       3.78        3.66
 Pearson's cor(R)  -0.74      0.53        -0.03
 Clean data (autosomes):
 %variance         21.4       3.82        3.66
 Pearson's cor(R)  -0.74      0.53        -0.07
 Best data (autosomes):
 %variance         22.6       3.7         3.6
 Pearson's cor(R)  -0.73      -0.25       -0.27