标签: biopython

我正在尝试更改以前的脚本，该脚本利用 biopython 获取有关物种门的信息。编写此脚本是为了一次检索一个物种的信息。我想修改脚本，以便我可以一次对 100 个生物执行此操作。这是初始代码

import sys
from Bio import Entrez

def get_tax_id(species):
    """to get data from ncbi taxomomy, we need to have the taxid.  we can
    get that by passing the species name to esearch, which will return
    the tax id"""
    species = species.replace(" ", "+").strip()
    search = Entrez.esearch(term = species, db = "taxonomy", retmode = "xml")
    record = Entrez.read(search)
    return record['IdList'][0]

def get_tax_data(taxid):
    """once we have the taxid, we can fetch the record"""
    search = Entrez.efetch(id = taxid, …

使用biopython库，我想删除列表中列出的残留物，如下所示。该线程（http://pelican.rsvs.ulaval.ca/mediawiki/index.php/Manipulated_PDB_files_using_BioPython）提供了一个去除残留物的示例。我有以下代码来去除残留物

 residue_ids_to_remove = [105, 5, 8, 10, 25, 48]
 structure = pdbparser.get_structure("3chy", "./3chy.pdb")
 first_model = structure[0]
 for chain in first_model:
     for residue in chain:
         id = residue.id
         if id[1] in residue_ids_to_remove:
             chain.detach_child(id[1])
 modified_first_model = first_model 

但这段代码不起作用并引发了错误

def detach_child(self, id):
    "Remove a child."
    child=self.child_dict[id]
    KeyError: '105'

这段代码有什么问题？

或者，我可以使用accept_residue()并将其写入PDB。我不想这样跟踪，因为我想在内存中执行此操作以进行进一步处理。

我正在编写一个 python 脚本，并且希望将查询序列信息作为字符串变量而不是 FASTA 格式文件传递到blastn（如果可能）。

我使用 Biopython 的 SeqIO 将多个转录本名称存储为键，并将其序列存储为关联值。

所以它看起来像这样

transcripts = dict()
for record in SeqIO.parse("transcript_sequences.fasta", "fasta"):
transcripts[record.name] = record.seq

print transcripts

所以字典看起来像这样

{'var_F': Seq('CTTCATTCTCGTTTAGCGGCTGCTCGTGGAAATTTCGAAAAAATCTGAAACTAG...TGC', SingleLetterAlphabet())}

现在我想将字典中的序列信息解析为爆炸查询和主题。

subprocess.call("blastn -query " + transcript['var_F'] + ".txt" + " -subject " + transcript['var_B'] + " -outfmt 6 > tmp_blast.txt", shell=True)  

我知道blast只接受文件名或字符串作为文件位置，但我只是想知道是否有解决方法。

预先感谢您阅读我的第一个问题：P

我有一些绳子，

['SGALWDV', 'GALWDVP', 'ALWDVPS', 'LWDVPSP', 'WDVPSPV']

这些字符串彼此部分重叠。如果您手动重叠它们，您将得到：

SGALWDVPSPV

我想要一种从重叠字符串列表到 python 中的最终压缩字符串的方法。我觉得这一定是一个已经有人解决的问题，并且正在努力避免重新发明轮子。我现在可以想象的方法要么是蛮力，要么是通过使用biopython和序列比对器变得比我想要的更复杂。我有一些简单的短字符串，只想以简单的方式正确合并它们。

有没有人有关于在 python 中执行此操作的好方法的建议？谢谢！

我有一个PDB文件'1abz'(https://files.rcsb.org/view/1ABZ.pdb),包含23个不同模型(编号为模型1-23)的蛋白质结构的坐标.请忽略标题备注,有趣的信息从第276行开始,标题为"模型1".

我想计算蛋白质的平均结构.蛋白质的PDB文件包含多个构象/模型,我想计算每个残基的单个原子的平均坐标,这样我最终得到一个构象/模型.

我无法弄清楚如何使用Biopython来做到这一点,所以我尝试使用Pandas计算平均坐标.我想我已经设法计算平均值,但现在的问题是我有一个不再是PDB格式的csv文件,所以我无法将此文件加载到PyMol中.

我的问题是,如何将我的csv文件转换为PDB格式.更好的是,如何在不损害原始pdb文件格式的情况下获取Biopython或Python中的平均坐标？

这是我用来计算熊猫平均坐标的代码.

#I first converted my pdb file to a csv file

import pandas as pd
import re

pdbfile = '1abz.pdb'
df = pd.DataFrame(columns=['Model','Residue','Seq','Atom','x','y','z']) #make dataframe object
i = 0 #counter

b = re.compile("MODEL\s+(\d+)")
regex1 = "([A-Z]+)\s+(\d+)\s+([^\s]+)\s+([A-Z]+)[+-]?\s+([A-Z]|)"
regex2 = "\s+(\d+)\s+([+-]?\d+\.\d+\s+[+-]?\d+\.\d+\s+[+-]?\d+\.\d+)"
reg = re.compile(regex1+regex2)

with open(pdbfile) as f:
    columns = ('label', 'ident', 'atomName', 'residue', 'chain', 'sequence', 'x', 'y', 'z', 'occ', 'temp', 'element')
    data = []
    for line in f:
        n = b.match(line)
        if n:
            modelNum = …

我是一名新近自学成才的（减去 1 节非常基础的课程）程序员，在生物实验室工作。我有一个脚本，它遍历来自两种不同细胞类型的 RNAseq 数据，并在另一个数据集中运行 ttest。它适用于这个应用程序，但代码感觉非常粗鲁，我知道我会写很多类似的脚本。

如何更好地编写以下代码以使其更高效？

计划目标：

1. 将基因列表与两种细胞类型的 rnaseq 文库进行比较，如果该文库包含该基因，则运行细胞类型 1 与细胞类型 2 的 ttest
2. 输出结果。

：

import pandas as pd
from scipy.stats import ttest_ind
rnatest = {'Gene symbol':["GeneA","GeneB"],"rnaseq1A":[1,1.5],"rnaseq1B":[1.3,1.2],"rnaseq2A":[2.3,2.7],"rnaseq2B":[2,2.6]} 
df = pd.DataFrame(rnatest)
GOIlist = ["GeneA","GeneB"]
GOI = []
mu = [] 
pval = []
for index, row in df.iterrows():
  if row['Gene symbol'] in GOIlist:
    t, p = ttest_ind([row["rnaseq1A"],row["rnaseq1B"]],[row["rnaseq2A"],row["rnaseq2B"]])
    GOI.append(row['Gene symbol'])
    mu.append(t)
    pval.append(p)
df2 = {'Gene symbol':GOI,"tVAL":mu, "pVAL":pval}
df2 = pd.DataFrame(df2)
print(df2)  

我有一个包含数千个登录号的文件：

看起来像这样..

>NC_033829.1 Kallithea virus isolate DrosEU46_Kharkiv_2014, complete genome
AGTCAGCAACGTCGATGTGGCGTACAATTTCTTGATTACATTTTTGTTCCTAACAAAATGTTGATATACT

>NC_020414.2 Escherichia phage UAB_Phi78, complete genome
TAGGCGTGTGTCAGGTCTCTCGGCCTCGGCCTCGCCGGGATGTCCCCATAGGGTGCCTGTGGGCGCTAGG

如果想将其拆分为多个文件，每个文件都有一个登录号，那么我可以使用以下代码

awk -F '|' '/^>/ {F=sprintf("%s.fasta",$2); print > F;next;} {print >> F;}' < yourfile.fa

我有一个包含数千个登录号（又名 >NC_*）的文件，并希望将其拆分，例如每个文件包含约 5000 个登录号。因为我是 awk/bash/python 的新手，所以我很难找到一个巧妙的解决方案

任何想法或评论表示赞赏

我有像这个小例子的文本文件:

>ENST00000491024.1|ENSG00000187583.6|OTTHUMG00000040756.4|OTTHUMT00000097942.2|PLEKHN1-003|PLEKHN1|176
SLESSPDAPDHTSETSHSPLYADPYTPPATSHRRVTDVRGLEEFLSAMQSARGPTPSSPLPSVPVSVPASDPRSCSSGPAGPYLLSKKGALQSRAAQRHRGSAKDGGPQPPDAPQLVSSAREGSPEPWLPLTDGRSPRRSRDPGYDHLWDETLSSSHQKCPQLGGPEASGGLVQWI
>ENST00000433179.2|ENSG00000187642.5|OTTHUMG00000040757.3|-|C1orf170-201|C1orf170|696
MPTQDGQLRRPARPPGPRAWMEPRGGGSSQFSSCPGPASSGDQMQRLLQGPAPRPPGEPPGSPKSPGHSTGSQRPPDSPGAPPRSPSRKKRRAVGAKGGGHTGASASAQTGSPLLPAASPETAKLMAKAGQEELGPGPAGAPEPGPRSPVQEDRPGPGLGLSTPVPVTEQGTDQIRTPRRAKLHTVSTTVWEALPDVSRAKSDMAVSTPASEPQPDRDMAVSTPASEPQSDRDMAVSTPASEPQPDTDMAVSTPASEPQPDRDMAVSIPASKPQSDTAVSTPASEPQSSVALSTPISKPQLDTDVAVSTPASKHGLDVALPTAGPVAKLEVASSPPVSEAVPRMTESSGLVSTPVPRADAAGLAWPPTRRAGPDVVEMEAVVSEPSAGAPGCCSGAPALGLTQVPRKKKVRFSVAGPSPNKPGSGQASARPSAPQTATGAHGGPGAWEAVAVGPRPHQPRILKHLPRPPPSAVTRVGPGSSFAVTLPEAYEFFFCDTIEENEEAEAAAAGQDPAGVQWPDMCEFFFPDVGAQRSRRRGSPEPLPRADPVPAPIPGDPVPISIPEVYEHFFFGEDRLEGVLGPAVPLPLQALEPPRSASEGAGPGTPLKPAVVERLHLALRRAGELRGPVPSFAFSQNDMCLVFVAFATWAVRTSDPHTPDAWKTALLANVGTISAIRYFRRQVGQGRRSHSPSPSS
>ENST00000341290.2|ENSG00000187642.5|OTTHUMG00000040757.3|OTTHUMT00000097943.2|C1orf170-001|C1orf170|676
MEPRGGGSSQFSSCPGPASSGDQMQRLLQGPAPRPPGEPPGSPKSPGHSTGSQRPPDSPGAPPRSPSRKKRRAVGAKGGGHTGASASAQTGSPLLPAASPETAKLMAKAGQEELGPGPAGAPEPGPRSPVQEDRPGPGLGLSTPVPVTEQGTDQIRTPRRAKLHTVSTTVWEALPDVSRAKSDMAVSTPASEPQPDRDMAVSTPASEPQSDRDMAVSTPASEPQPDTDMAVSTPASEPQPDRDMAVSIPASKPQSDTAVSTPASEPQSSVALSTPISKPQLDTDVAVSTPASKHGLDVALPTAGPVAKLEVASSPPVSEAVPRMTESSGLVSTPVPRADAAGLAWPPTRRAGPDVVEMEAVVSEPSAGAPGCCSGAPALGLTQVPRKKKVRFSVAGPSPNKPGSGQASARPSAPQTATGAHGGPGAWEAVAVGPRPHQPRILKHLPRPPPSAVTRVGPGSSFAVTLPEAYEFFFCDTIEENEEAEAAAAGQDPAGVQWPDMCEFFFPDVGAQRSRRRGSPEPLPRADPVPAPIPGDPVPISIPEVYEHFFFGEDRLEGVLGPAVPLPLQALEPPRSASEGAGPGTPLKPAVVERLHLALRRAGELRGPVPSFAFSQNDMCLVFVAFATWAVRTSDPHTPDAWKTALLANVGTISAIRYFRRQVGQGRRSHSPSPSS
>ENST00000428771.2|ENSG00000188290.6|OTTHUMG00000040758.2|OTTHUMT00000097945.2|HES4-002|HES4|247
MAADTPGKPSASPMAGAPASASRTPDKPRSAAEHRKVGSRPGVRGATGGREGRGTQPVPDPQSSKPVMEKRRRARINESLAQLKTLILDALRKESSRHSKLEKADILEMTVRHLRSLRRVQVTAALSADPAVLGKYRAGFHECLAEVNRFLAGCEGVPADVRSRLLGHLAACLRQLGPSRRPASLSPAAPAEAPAPEVYAGRPLLPSLGGPFPLLAPPLLPGLTRALPAAPRAGPQGPGGPWRPWLR

该文件被分成不同的组.每组有2个部分.第1部分开始,">"并且该部分中的元素被分割,"|"并且之后的线是第2部分.我试图从我的文件中创建一个Python列表,其中包含每个组ID部分的第6个元素.以下是小例子的预期输出:

list = ["PLEKHN1", "C1orf170", "C1orf170", "HES4"]

我试图先导入一个字典,然后使用以下方法创建一个像预期输出的列表:

from itertools import groupby
with open('infile.txt') as f:
    groups = groupby(f, key=lambda x: not x.startswith(">"))
    d = {}
    for k,v in groups:
        if not k:
            key, val = list(v)[0].rstrip(), "".join(map(str.rstrip,next(groups)[1],""))
            d[key] = val

k = d.keys()
res = [el[5:] for s in k for el in s.split('|')]

但它不会返回我想要的东西.你知道怎么解决吗？

1. 有没有办法将 hmmer 输出转换为 pandas 数据帧？
2. 我也不确定如何通过 Bio 模块将 hmmer tblout 表加载到 python 中。

我相信您可以使用 SeqIO.parse 或 SeqIO.search 调用 hmmer 格式。表格的格式显示为制表符分隔，但它似乎是随机空格的集合，这意味着如果我删除标题并 # 只留下表格信息使用制表符分隔符分割表格并不简单。

hmmer --tblout 文件的一个小示例如下：

#                                                                                       --- full sequence ---- --- best 1 domain ---- --- domain number estimation ----
# target name        accession  query name                                   accession    E-value  score  bias   E-value  score  bias   exp reg clu  ov env dom rep inc description of target
#------------------- ----------                         -------------------- ---------- --------- ------ ----- --------- ------ -----   --- --- --- --- --- --- …

如果之前有人问过这个问题，我深表歉意，但我已经搜索了好几天，但在 Python 中找不到解决方案。

我有一个很大的 fasta 文件，包含标题和序列。

>cavPor3_rmsk_tRNA-Leu-TTA(m) range=chrM:2643-2717 5'pad=0 3'pad=0 strand=+ repeatMasking=none
GTTAAGGTGGCAGAGCCGGTAATTGCATAAAATTTAAGACTTTACTCTCA
GAGGTTCAACTCCTCTCCTTAACAC

>cavPor3_rmsk_tRNA-Gln-CAA_ range=chrM:3745-3815 5'pad=0 3'pad=0 strand=- repeatMasking=none
AGAGGGTCATAAAGGTTATGGGGTTGGCTTGAAACCAGCTTTAGGGGGTT
CAATTCCTTCCTCTCT

>cavPor3_rmsk_tRNA-Ser-TCA(m) range=chrM:6875-6940 5'pad=0 3'pad=0 strand=- repeatMasking=none
AGAGGGTCATAAAGGTTATGGGGTTGGCTTGAAACCAGCTTTAGGGGGTT
CAATTCCTTCCTCTCT

这是文件的一个非常小的片段。如果正如您在最后两个条目中看到的那样，序列是相同的，我只想保留第一个条目（标题和序列）。

输出如下所示：

>cavPor3_rmsk_tRNA-Leu-TTA(m) range=chrM:2643-2717 5'pad=0 3'pad=0 strand=+ repeatMasking=none
GTTAAGGTGGCAGAGCCGGTAATTGCATAAAATTTAAGACTTTACTCTCA
GAGGTTCAACTCCTCTCCTTAACAC

>cavPor3_rmsk_tRNA-Gln-CAA_ range=chrM:3745-3815 5'pad=0 3'pad=0 strand=- repeatMasking=none
AGAGGGTCATAAAGGTTATGGGGTTGGCTTGAAACCAGCTTTAGGGGGTT
CAATTCCTTCCTCTCT

问题是 FASTA 文件大小超过 1 GB。我已经找到了解决这个问题的方法，方法是根据重复的 ID 删除重复项或使用 bash，但遗憾的是我无法在我的计算机上执行此操作。此任务是一个研究项目，而不是家庭作业或任务。

预先感谢您的帮助！