!head /Volumes/web/cnidarian/TJGR_v9_cpgreport_1.output

!head /Volumes/web/cnidarian/TJGR_v9_cpgreport_1.gff

!head /Volumes/web/cnidarian/TJGR_v9_newcpgreport_1.output

Failed
./cpgplot -sequence /Volumes/Bay3/Software/bismark_v0.6.4/genomes/oyster.v9.fa

!head /Volumes/web/cnidarian/TJGR_v9cpgplot1.gff


!fgrep -c "sequence_feature" /Volumes/web/cnidarian/TJGR_v9cpgplot1.gff

!head /Volumes/web/cnidarian/TJGR_v9cpgplot2.gff


!fgrep -c "sequence_feature" /Volumes/web/cnidarian/TJGR_v9cpgplot2.gff

!head /Volumes/web/cnidarian/TJGR_v9cpgplot3.gff


!fgrep -c "sequence_feature" /Volumes/web/cnidarian/TJGR_v9cpgplot3.gff

#how about intersecting methylation with CpG islands and plot histogram?
!intersectbed -c -a /Volumes/web/cnidarian/TJGR_v9cpgplot3.gff -b /Volumes/web/bivalvia/wholegenomefiles_MBDbsSeq_gill/finalSupp_figshare/Cgigas_gill_HTbisulfiteSeq_CG_methylated.gff > /Volumes/web/cnidarian/BiGill_island3_intersect_mCpG.txt
!head /Volumes/web/cnidarian/BiGill_island3_intersect_mCpG.txt

!wc /Volumes/web/cnidarian/BiGill_island3_intersect_mCpG.txt

from pandas import *

# read data from data file into a pandas DataFrame  
island = read_table("/Volumes/web/cnidarian/BiGill_island3_intersect_mCpG.txt", # name of the data file
            #sep="\t", # what character separates each column?
            #na_values=["", " "], # what values should be considered "blank" values?
            header=None)



island[9].hist(bins=50);
#Axis limits are changed using the axis([xmin, xmax, ymin, ymax]) function.
#plt.axis([0, 30, 0, 1000]);
plt.title('island_meth');

#also need to get total CG...

!intersectbed -c -a /Volumes/web/cnidarian/TJGR_v9cpgplot3.gff -b /Volumes/web/trilobite/Crassostrea_gigas_v9_tracks/Cgigas_v9_CG.gff > /Volumes/web/cnidarian/BiGill_island3_intersect_CG.txt

!head /Volumes/web/cnidarian/BiGill_island3_intersect_CG.txt

!wc /Volumes/web/cnidarian/BiGill_island3_intersect_CG.txt

# read data from data file into a pandas DataFrame  
iscg = read_table("/Volumes/web/cnidarian/BiGill_island3_intersect_CG.txt", # name of the data file
            #sep="\t", # what character separates each column?
            #na_values=["", " "], # what values should be considered "blank" values?
            header=None)

iscg[9].hist(bins=50);
#Axis limits are changed using the axis([xmin, xmax, ymin, ymax]) function.
#plt.axis([0, 1000, 0, 200]);
plt.title('island_cg');

#now need to do some math (SQLShare?) for each island mCG/totalCG to get %methylation

print island[9],iscg[9],island[9]/iscg[9]



GAVE UP -- excelled it


!head /Volumes/web/cnidarian/BiGill_island3_intersect_PERmCpG.txt

# read data from data file into a pandas DataFrame  
perMcg = read_table("/Volumes/web/cnidarian/BiGill_island3_intersect_PERmCpG.txt", # name of the data file
            #sep="\t", # what character separates each column?
            #na_values=["", " "], # what values should be considered "blank" values?
            header=None)

perMcg[11].hist(bins=50);
#Axis limits are changed using the axis([xmin, xmax, ymin, ymax]) function.
plt.axis([0, 1, 0, 2000]);
plt.title('Island percent_mcg');

# Develop track of CpG islands that are >80% methylated

!head /Volumes/web/cnidarian/TJGR_v9cpgplot3_80.gff

!wc /Volumes/web/cnidarian/TJGR_v9cpgplot3_80.gff

#ID promoters that do not overlap gene bodies



!subtractbed -a /Volumes/web/trilobite/Crassostrea_gigas_v9_tracks/Cgigas_v9_1k5p_gene_promoter.gff -b /Volumes/web/trilobite/Crassostrea_gigas_v9_tracks/Cgigas_v9_gene.gff > /Volumes/web/cnidarian/TJGR_prom_subtract_gene1.gff 

!head /Volumes/web/cnidarian/TJGR_prom_subtract_gene1.gff 

!tail /Volumes/web/cnidarian/TJGR_prom_subtract_gene1.gff 

!wc /Volumes/web/cnidarian/TJGR_prom_subtract_gene1.gff 

!wc /Volumes/web/trilobite/Crassostrea_gigas_v9_tracks/Cgigas_v9_1k5p_gene_promoter.gff



!intersectbed -wa -a /Volumes/web/cnidarian/TJGR_prom_subtract_gene1.gff -b /Volumes/web/cnidarian/TJGR_v9cpgplot3.gff > /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1.gff  

!head /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1.gff  

!tail /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1.gff  
#remove duplicate lines (occurred when more than one island present)
!uniq  /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1.gff  > /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1u.gff
!wc /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1.gff  

!wc /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1u.gff  

!tail /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1u.gff  

#again with more generous definition of CpG Island 

!intersectbed -wa -a /Volumes/web/cnidarian/TJGR_prom_subtract_gene1.gff -b /Volumes/web/cnidarian/TJGR_v9cpgplot2.gff > /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2.gff  

#removing duplicates lines
!uniq  /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2.gff  > /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u.gff

!wc /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2.gff  

!wc /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u.gff  

!head /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u.gff 

!intersectbed -c -a /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u.gff  -b /Volumes/web/trilobite/Crassostrea_gigas_v9_tracks/Cgigas_v9_CG.gff > /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG.gff 

!head /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG.gff

!intersectbed -a /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u.gff  -b /Volumes/web/trilobite/Crassostrea_gigas_v9_tracks/Cgigas_v9_CG.gff -c > /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG2.gff  

!tail /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG2.gff  

!sed 's/Roberts_20100712_CC_F3_trimmed/Haliotis_cra_v3/g' </Volumes/web/cnidarian/lft_BlackAbalone_v3_swissprot_blastout_c> /Volumes/web/cnidarian/lft_BlackAbalone_v3_swissprot_blastout_d
#sed 's/abc/XYZ/g' <infile >outfile
!tr "\n \t" ' ' </Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG2.gff> /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG3.gff

!head /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland2u_in_CG3.gff

#absolute methylation
#using the conserved model
#Window size [100]: 
#Minimum length of an island [200]: 
#Minimum observed/expected [0.6]: 
#Minimum percentage [50.]: 45
#promoter does not overlap gene
!intersectbed -c -a /Volumes/web/cnidarian/TJGR_prom_notgene_cpgIsland1u.gff -b /Volumes/web/bivalvia/wholegenomefiles_MBDbsSeq_gill/finalSupp_figshare/Cgigas_gill_HTbisulfiteSeq_CG_methylated.gff > /Volumes/web/cnidarian/BiGill_pro_island1u_intersect_mCpG.txt

!head /Volumes/web/cnidarian/BiGill_pro_island1u_intersect_mCpG_mod.txt

!wc /Volumes/web/cnidarian/BiGill_pro_island1u_intersect_mCpG_mod.txt

from pandas import *
# read data from data file into a pandas DataFrame  
npmpg = read_table("/Volumes/web/cnidarian/BiGill_pro_island1u_intersect_mCpG_mod.txt", # name of the data file
            #sep="\t", # what character separates each column?
            #na_values=["", " "], # what values should be considered "blank" values?
            header=None)

npmpg[9].hist(bins=50);
#Axis limits are changed using the axis([xmin, xmax, ymin, ymax]) function.
#plt.axis([0, 40, 0, 500]);
plt.title('nProm_mCpG');

#next step plot absolute methylation against gene expression as provided by Mac