Noi/NOTES/2012-7-28

From ZhangLabWiki
Jump to navigation Jump to search

[Link to calendar]

overlapping samples for Illumina 450k and BSPP (GA data set)[edit]

  • Aim:
    • To get the CpG sites that overlap between the two methods, Illumina 450k methylation arrays and BSPP assay, then calculate the correlation coefficient (will be performed by Greg)
    • To extract DNA methylation values of the CpG sites in the genes associated with the 69 predictive markets (as shown in Table_S3).

35 sample lists[edit]

SampleID MethylFreq Corr
(minReadDepth=10)
MethylFreq Corr
(minReadDepth=20)
773002 0.956 0.967
1069001 0.961 0.970
1100001 0.959 0.970
825002 0.939 0.955
880003 0.935 0.952
888003 0.938 0.956
1373001 0.932 0.951
6925001 0.932 0.953
6988001 0.845 0.880
7161001 0.932 0.952
7419001 0.883 0.915
7528001 0.966 0.976
7529001 0.936 0.954
7729002 0.887 0.917
7781002 0.752 0.761
10174001 0.890 0.922
10182001 0.828 0.843
10272001 0.962 0.972
10703001 0.851 0.870
10741001 0.893 0.921
10914001 0.864 0.903
11029001 0.849 0.884
11066001 0.919 0.942
11200001 0.902 0.928
11202001 0.885 0.917
11285001 0.932 0.952
11324001 0.969 0.977
11340001 0.919 0.940
11379001 0.929 0.950
11394002 0.938 0.956
11434001 0.891 0.922
11594001 0.927 0.945
11595001 0.961 0.970
11756001 0.970 0.977
11793001 0.966 0.974
In genome-miner
Directory: /home/nplongth/Noi_scratch/BSPP-AMD_20111012/AMD_mappig.data_hg19/AMD-35samples
BED files of 35 samples are in sub directory -> GA35samples
Illumina 450k marker is in Illumina450kMarkers.txt
# I generate methylFreq matrix (minRD 10, min_sample 5 and min_sample28 (80%))
## minRD 10, min_sample 5: total number of CpG sites: 693,202 sites (GA-35samplemethyl_min5_minSTD0)
## minRD 10, min_sample 28: total number of CpG sites: 413,904 sites (GA-35samplemethyl_RD10_min28_minSTD0)

# minRD 10, min_sample 5, number of sites overlapping with 450k : 61,296 sites (Illumina-450k-overlapped_GA-35samples-minRD10_minSample1.txt)
# minRD 10, min_sample 28, number of sites overlapping with 450k : 37,239 sites (Illumina-450k-overlapped_GA-35samples-minRD10_minSample28.txt)

These numbers showed very low overlapping of the CpG sites covered by the two assays. I asked Greg to confirm that the coordinates are reported to hg19.
I also check when the coordinate is -1, -2 or +1 or -2 (change the coordinate of Illuminar markers (chr_pos),
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 <Illu-1.txt | wc -l -> 0
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 <Illu-2.txt | wc -l -> 3,112 ?
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 < Illu+1.txt | wc -l -> 0
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 < Illu+2.txt | wc -l -> 3,104 ?
This is unlikely to get more overlapping from both -2 and +2 . However, I still need to check carefully.
# check duplication among the site +2, -2 to the overlapping sites from original coordinate.
awk '{print $1}' Illumina-450k-overlapped_GA-35samples-minRD10_minSample5.txt > tmp-ori (61,295 sites)
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 < Illu-2.txt | awk '{print $1}' > tmp-2 (3,112 sites)
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 < Illu+2.txt | awk '{print $1}' > tmp+2 (3,104 sies)
cat tmp-ori tmp-2 | sort | uniq -d | wc -l -> 3,009 
cat tmp-ori tmp-2 | sort | uniq -D | wc -l -> 6,018 

cat tmp-ori tmp+2 | sort | uniq -d | wc -l ->3,016
cat tmp-ori tmp+2 | sort | uniq -D | wc -l -> 6,032
This suggests that the sites that I saw when shift -2 are just the sites that have very close position and already listed in the overlapping sites. I will ignore ~ 100 sites from each +2 and -2 since it is unlikely to be to right overlapping sites. 
I then look at UCLA data to see if the CpG sites covered are consistently covered in GA samples and overlapping with 450k arrays markers.
# I generated matrix (minRD 10, min_sample 5: UCLA_minsample5_minSTD0 -> 691,556 sites total
# number of sites overlapping with 450k: 61,273 sites
# number of sites overlapping with GA data set: 649,623 (94%)

Discussion with Dr. Zhang[edit]

  • 2012.07.30

Extract the CpG position in the genes associated with the 69 predictive markers[edit]

  1. Download Human RefSeq from genomebrowser: http://genome.ucsc.edu/cgi-bin/hgTables -> output file: Hg19.Refseq.gz -> uncompressed -> Hg19.Refseq

Header of Hg19.Refseq:
#bin(1) name(2) chrom(3) strand(4) txStart(5) txEnd(6) cdsStart(7) cdsEnd(8) exonCount(9) exonStarts(10) exonEnds(11) score(12) name2(13) cdsStartStat(14) cdsEndStat(15) exonFrames(16)

 head -10 Hg19.Refseq | sed 's/\t/\n/g' 
header:
(1)#bin
(2)name
(3)chrom
(4)strand
(5)txStart
(6)txEnd
(7)cdsStart
(8)cdsEnd
(9)exonCount
(10)exonStarts
(11)exonEnds
(12)score
(13)name2
(14)cdsStartStat
(15)cdsEndStat
(16)exonFrames

To generate bed file containing Chromosome position, strart position of the gene (-2,000 bp window), end postion of genes (+2,000 bp window), and gene names
awk '{print $3"\t"$5-2000"\t"$6+2000"\t"$13}' Hg19.Refseq | sort -k1,1 -k2,2n -u > Hg19.Refseq.bed

To fix the negative integers of start position
awk '{if($2 < 0 ) $2=0; print $0}' Hg19.Refseq.bed > Hg19.Refseq.bed.fixed
/home/dinh/softwares/BEDTools-Version-2.16.2/bin/bedtools intersect -wao -a GA-35samples_minRD10_min5.BED.txt -b Hg19.Refseq.bed.fixed | awk '{if ($9 ~ /1/) print $0}' > GA-35samples_minRD10_min5-Hg19Refseq2kbwindow.fixed
./get_overlapped_site.pl ori.tableS3.single-genelist.txt GA-35samples_minRD10_min5-Hg19Refseq2kbwindow.fixed > GA-35samples_minRD10_min5-Hg19Refseq2kbwindow-69markers-OL
./extract-35samples-69Markers.pl GA-35samplemethyl_min5_minSTD0 <GA-35samples_minRD10_min5-Hg19Refseq2kbwindow-69markers-OL| sort |uniq -u > methylData.minRD10_minsmaple5.2kbwindow.69markers-assocgenes.txt
./extract-35samplesCpG.pl GA-35samplemethyl_min5_minSTD0 < tableS3-2.txt > methylData.minRD10_minsmaple5.69markers-overlap.txt

Rename the files and send to Greg

  1. GA-35samples_minRD10_minsample5-2kbwindow-69markers-assocgenes.gz
  2. methylData.minRD10_minsmaple5.2kbwindow.69markers-assocgenes.txt.gz
  3. methylData.minRD10_minsmaple5.69markers-overlap.txt.gz