Calculate cumulative gc skew in python

Question

[bird77]

I want to write a python function that calculates the cumulative gc skew along a entire bacterial genome.

I have no problem implementing GC skew (G-C)/(G+C) but what would I need to change to calculate the cumulative GC skew?

Thank you.

[Joe]

You want to do something like calculate the GC skew across a window of desired size, then slide that window accross the genome? Why cumulative, out of interest?

Do you need to write this from scratch (e.g. for homework), or are tools such as BioPython's implementation an option?

[bird77]

I thought this is the method for finding the origin of replication?

Answer 1

[Joe]

Maybe this will give you something to work with...

Given a fasta input sequence,

from Bio import SeqIO, SeqUtils

rec = SeqIO.read('/path/to/seqs.fasta', 'fasta')
gc = SeqUtils.GC_skew(rec.seq, 100)

>>> gc
[-0.2682926829268293, -0.05660377358490566, -0.6363636363636364, -0.2, -0.20689655172413793, -0.2682926829268293, 0.01818181818181818, -0.7222222222222222, -0.3333333333333333, -0.44, -0.2558139534883721, -0.37254901960784315, -0.35714285714285715, -0.3877551020408163, -0.52, -0.44, -0.4583333333333333, -0.44680851063829785, -0.6666666666666666, -0.5636363636363636, -0.24, -0.07692307692307693, -0.19148936170212766, -0.12280701754385964, -0.30434782608695654, -0.5121951219512195, -0.23636363636363636, -0.48, 0.0, -0.10204081632653061, -0.09433962264150944, -0.19047619047619047, -0.40540540540540543, -0.2413793103448276, -0.13636363636363635, -0.37254901960784315, -0.25, -0.2571428571428571, -0.375, -0.3333333333333333, -0.20689655172413793, -0.19298245614035087, -0.3488372093023256, -0.36363636363636365, -0.14285714285714285]

# There's also a graphical output

SeqUtils.xGC_skew(rec.seq, 100)

Which produces (not that it's the most useful plot in the world):

To use the accumulated GC for your own purposes, just sum the gc list.