I have the assembly file of a certain lizard species, I wanted to recreate Fig. 2 of the following paper, which compares the GC content distribution over 5kb windows, how would I go about doing this?
https://bmcecolevol.biomedcentral.com/articles/10.1186/s12862-022-02086-7#Sec10
I would suggest using bedtools makewindows and bedtools nuc to compute the GC content for genomic windows as shown in Fig 2 of the linked paper.
See related question here.
This will slice the genome into 5Kbp (or whatever size specified to makewindows) pieces and compute the GC content of each piece, which can be plotted as a distribution as in the figure.
Something like this would work for the 5000 bp example:
# assemble a bedtools compatible genome file:
# see https://www.biostars.org/p/70795/ for hints how to do this.
bedtools makewindows -w 5000 -g my_genome.genome > windows.bed # you just made this .genome file manually
bedtools nuc -bed windows.bed -fi my_genome.fasta > windows_gc.bed
The bed output file will contain the GC content of each window in the 7_pct_gc column (according to my understanding of the output file).
samtools. You could simply generate an index file by running samtools faidx genome.fasta, and use the output file with extension .fasta.fai as input in your makewindwos command.
$\endgroup$
Commented
Sep 16, 2023 at 13:34
tail -n +2 windows_gc.bed | cut -f 5 | awk '{printf("%d\n",$1*100)}' | sort | uniq -c | sed 's/^ *//'| sed 's/ /,/'>output.csv. This should print a CSV file which you can easily plot.
$\endgroup$
Commented
Sep 16, 2023 at 13:59
Try this small package/library: https://github.com/willros/GuCi/
It produces a plot of the gc content of the whole genome. It works on small genomes, like bacteria.
Hope you like!
