SNP density per genes
0
0
Entering edit mode
2.3 years ago
rthapa ▴ 50

Hi,

I want to estimate SNP density per gene (>20,000 genes). I do have a gff file and a vcf file with SNPs position. I am trying to use the following python script but couldn't get the density. I am new to python so there could be something wrong with my script. Does anyone have any suggestions?

Thanks

import allel
allel.__version__

geneset = allel.gff3_to_dataframe('gene.gff3', attributes=['Name'])
geneset.head()

genes = geneset[geneset['type'] == 'gene']

callset = allel.read_vcf("SAP.vcf")
gt = allel.GenotypeArray(callset['calldata/GT'])

import numpy as np
import numcodecs
import sys
n_genes = len(genes)

for i, (_, gene) in genes.iterrows():
    if geneset[0] == 'Chr01'
        try:
            loc_gene = pos.locate_range(gene.start, gene.end)
        except KeyError:
            pass
        else:
            n_variants[i] = loc_gene.stop - loc_gene.start
SNP • 396 views
ADD COMMENT

Login before adding your answer.

Traffic: 3110 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6