Generate distribution of conserved noncoding elements and noncoding spacers of Drosophila genomes

[MManee]

• we need to identify:
  • all coding regions and all noncoding regions
  • conserved elements and spacers in coding regions
  • conserved elements and spacers in nonoding regions
• lets start with Drosophila melanogaster
  • we will need three datasets:
  • most conserved elements in Drosophila genomes (phastCons), and can be downloaded from http://hgdownload.soe.ucsc.edu/goldenPath/dm3/database/phastConsElements15way.txt.gz
  • exons in D. melanogaster genome, and it is located here: https://github.com/kacst-bioinfo-lab/labwork/tree/master/data/shared
  • chromosome coordinates of D. melanogaster, and can be downloaded from: http://hgdownload.soe.ucsc.edu/goldenPath/dm3/database/chromInfo.txt.gz

[MManee]

for newbies in GitHub, this thread #1 , should be useful!

[amer1404]

For coding steps :

• add starting column (0) for each gene to "chromInfo.txt"
• rename "chromInfo.txt" to "genes.bed"
• rename "dm3_exons.txt" to reads.bed"
• from bedtools, I've run the intersect command :

bedtools intersect -a reads.bed -b genes.bed > coding.bed For non coding steps:

• create new genes file "genes_non.bed" to fit bedtools format., the new genes file contain only name of gene and end of gene.
• run the complement command:

bedtools complement -i coding.bed -g genes_non.bed > non-coding.bed

• after running the complement command appear the following error:

  • Error: Sorted input specified, but the file coding.bed has the following out of order record

    chr4    252055    252253    CG1674-RC_exon_0_0_chr4_252056_f    0    +

  • to solve this problem, I've sorted the "coding.bed" records using bedtools command:

    sortBed -i coding.bed > sorted_coding.bed

  • as a result of sorting command we got new sorted file which is "sorted_coding.bed"

• then, I've running the complement command with new sorted file , but I got new error which is:

  • Error: Sorted input specified, but the file sorted_coding.bed has the following record with a different sort order than the genomeFile chromInfo.bed

    chrU    18133    18346    CG40189-RA_exon_0_0_chrU_18134_f    0    +

• to solve new error, I've merged the "sorted_coding.bed" records using bedtools comand bedtools merge sorted_coding.bed -i > merged_sorted_coding.bed

• then, I've running the complement command with new merge sorted file "merged_sorted_coding.bed" , but I got same previous erroe.

  bedtools complement -i merged_sorted_coding.bed -g genes_non.bed > non-coding.bed

• to solve this error, I tracking the error by return to "genes_non.bed" file. then, I move the "chrU" recored to end of file then save file and try to rerun the complement command again.

• after running the complement command, I got same error with new record which is "chrUextra" .

• I repeated the previous step with this erroe by moving "chrUextra" at the end of "genes_non.bed" file.

• also I got same error with next record , so I repeated the previous step with remaining records. after this step the error was going

• please find these files:

  • "coding.bed" for intersect exons.
  • "non-coding.bed" for complement exons.
  • "non-coding_merged.bed" for complement exons after merged the exons records

[MManee]

• genome browser can be used to validate the custom tracks you generated: https://genome.ucsc.edu/cgi-bin/hgCustom
• another way to validate your generated tracks is to perform the following calculations:
  • the total length of coding regions (bp) plus the total length of noncoding regions (bp) should be the total length of the genome (bp)
• use both of them to confirm the correctness of your results

[amer1404]

To make sure the results of previous files, I've calculate of coding lengths and noncoding lengths and compared with chromosome lengths. Results of calculate were not equal.

to solve this problem, We have to try the following steps:

1- find coding and non coding data from orginal files using bedtools.

coding_subtract.txt non_coding_subtract.txt bedtools intersect -a reads.bed -b genes.bed > coding.bed 2- find non coding data using bedtools. bedtools complement -i coding.bed -g genes_non.bed > non-coding.bed 3- we consider solution in previous issues that solve complement error.

4- find the coding data agin form noncoding data by run complement command 5- run subtract command on both files Rresult of calculation :

total non coding length :139034833 total coding length :29701678 total coding +non coding :168736511

gene length :168736537

[MManee]

• If you want to calculate the length of, for example, coding regions:

  genomeCoverageBed -i coding.bed -g allgenome 

[MManee]

• steps to generate coding regions an noncoding regions of Drosophila melanogaster:

  • sort chromosmes (all genome) and exons (coding regions) files using this:

    sort -k 1,1 -k2,2n allgenomefile > sorted_genome
    sort -k 1,1 -k2,2n dm3_exons.txt > sorted_dm3_exons.bed

  • use bedtools complement to generate noncoding regions

    bedtools complement -i sorted_dm3_exons.bed -g sorted_genome > noncoding.bed

  • the bp lengths of coding and noncoding regions were calculated:

  category	length (bp)
  coding regions	29701704
  noncoding regions	139034833
  total	168736537

  • coding and noncoding datasets can be found here: https://github.com/kacst-bioinfo-lab/labwork/tree/master/data/manee
  • next step is to split noncoding regions into intronic and intergenic regions

[amer1404]

please find the intronic data and perl code

code+data.tar.gz

[amer1404]

last update data.tar.gz

[amer1404]

this is final steps and result of Perl code : 1- downloded chromosom ( genom) ==> sorted_genome.txt 2- sort the allgenomefile.txt into sorted_genome.txt sort -k 1,1 -k2,2n allgenomefile.txt > sorted_genome.txt 3- sort the dm3_exons.txt to sorted_dm3_exons.bed sort -k 1,1 -k2,2n dm3_exons.txt > sorted_dm3_exons.bed 4- creat new genome file with valid entries by perl 5- creat the noncoding file bedtools complement -i sorted_dm3_exons.bed -g sorted_genome > noncoding.bed 6- creat the coding file by copying sorted_dm3_exons.bed to coding.bed 7- added 0 as start last genom file by perl code 8- bedtools subtract -a sorted_genome.bed -b noncoding.bed > coding.bed 9- create intronic file using cooding file without +1 -1 10- system("bedtools subtract -a ../data/noncoding.bed -b ../data/intronic.bed > ../data/intergenic.bed"); 11- calculate lengh of intronic + intergenic = noncoding lenght 12- calculate lengh of noncoding + codind = total lenght of genome

Results :

total lenght of coding                   : 29701704
total lenght of noncoding                : 139034833
 total of lenght of noncoding and coding : 168736537 
total lenght of genome                   : 168736537

total lenght of intronic                : 109319658
 total lenght of intergenic             : 29715175
 total lenght of intergenic and intronic: 139034833
 total lenght of noncoding              : 139034833

labwork.tar.gz

[MManee]

This is so awesome! Well done, Amer. The next step is to generate conserved noncoding elements (CNEs) and noncoding spacers of intronic and intergenic regions.

[atalgarni]

Alrighty, I just started ... so far I have just edited the chromInfo file. The Code can be found at https://github.com/kacst-bioinfo-lab/labwork/tree/master/code/abdulmalek/drosphila.py To make things easier and reproducible I left the input to the user as seen in lines 7-8,,however I am still stuck at how to automate this with os.system at lines 30-31 For the algorithm: 1- Add new column with 0 input for all rows inside chromInfo.txt, and convert the file into bed format. 2- Sort the chromInfo.bed using system sort | uniq functions. 3- Sort the dm3_exons.txt using sort | uniq functions.

For the output file please check https://github.com/kacst-bioinfo-lab/labwork/tree/master/doc/abdulmalek/chromInfo.bed https://github.com/kacst-bioinfo-lab/labwork/tree/master/doc/abdulmalek/chromSort.bed https://github.com/kacst-bioinfo-lab/labwork/tree/master/doc/abdulmalek/dm3_exonsSort.bed

@MManee @amer1404 @SULTAN-ALHARBI

Cheers,

[Maali055]

Print "Hello Group! \n";

I would like to thank Dr. Manee for involving me in the group. I Hope to be a positive addition to the group and absolutely I will learn alot from you.

[MManee]

• good related paper to read: Dynamic Epigenetic Control of Highly Conserved Noncoding Elements http://journals.plos.org/plosone/article?id=10.1371/journal.pone.0109326

[amer1404]

Hello Everyone

please find below the final result, Tomorrow I will upload all related files.

Total lenght of coding                         : 29701704
Total lenght of noncoding                      : 139034833
Total of lenght of noncoding and coding        : 168736537 
Total lenght of genome                         : 168736537

Total lenght of intronic                       : 109319658
Total lenght of intergenic                     : 29715175
Total lenght of intergenic and intronic        : 139034833
Total lenght of noncoding                      : 139034833

*************************************************************************************
Total lenght of coding CEs                    : 18219044
Total lenght of coding spacer                  : 11482662
Total lenght of CEs and spacer coding         : 29701706
Total lenght of coding                         : 29701704

*************************************************************************************
Total lenght of intergenic CNEs                : 2744358
Total lenght of intergenic spacer              : 26970817
Total lenght of CNEs and spacer coding         : 29715175
Total lenght of intergenic                     : 29715175

*************************************************************************************
Total lenght of intronic CNEs                  : 31527418
Total lenght of intronic spacer                : 77792240
Total lenght of CNEs and spacer coding         : 109319658
Total lenght of intronic                       : 109319658

*************************************************************************************`

[MManee]

It looks perfect. We will need to check the datasets using genome browser. Well done Amer.

[amer1404]

Hello dear All, :)

Kindly please find the updated files in my folders( data and code).