标签: bioinformatics

我有一个缩写列表

字母= ['Ala','Asx','Cys',......'Glx']

我想将其输出到一个文本文件,看起来像这样:

#Letters
Ala,Asx,Cys,..... Glx

Noob程序员在这里!我总是忘记最简单的事情!啊

请帮助和谢谢!

import Bio
from Bio import Seq
from Bio.Seq import Alphabet

output = 'alphabetSoupOutput.txt'
fh = open(output, 'w')
ThreeLetterProtein = '#Three Letter Protein'
Letters = Bio.Alphabet.ThreeLetterProtein.letters
fh.write(ThreeLetterProtein + '\n')

  #Don't know what goes here

fh.close()

我正在开发一个java程序,其中我正在使用一大类Class实例.我必须为数组的每个对象进行一些计算(使用它的Class方法).当我使用Arrays.Sort或其他方法对数组进行排序时,我注意到计算时间的增加以执行我的计算,即排除进行排序的时间!

如果我对Arrays.sort函数进行注释,则在其他计算之前,计算时间约为1或2毫秒,而取消注释则增加到5或6毫秒.

我必须为大量的步骤(甚至300,000)做这项工作,所以增量会改变很多我的表现.

对这种奇怪的行为有什么解释吗？

所以我有一个递归代码,可以为2条DNA链提供最佳的对齐,但问题是它的执行速度非常慢(我需要它递归).然后我在麻省理工学院网站上看到结果是附加的,这对我很好,但后来我想了一下,我发现有一个问题:

网站:http://ocw.mit.edu/courses/electrical-engineering-and-computer-science/6-096-algorithms-for-computational-biology-spring-2005/lecture-notes/lecture5_newest.pdf

麻省理工学院网站说,对于给定的溢出(i,j):first_strand(0,i)和second_strand(0,j)alignment
+ first_strand(i,len)和second_strand(j,len)对齐

等于

first_strand和second strand alignment

但:

GTC GTAA

具有GTA比对的G是G--并且具有A比对的GTA TC是TC并且A-结果= G-TC和GTAA-

真正的最佳结果= GTC-GTAA

谁能解释他们在麻省理工学院网站上的意思？我可能错了!

目前我正在努力应对我尚未能解决的AWK问题.我有一个巨大的文件(30GB),基因组数据包含一个包含位置的列表(在第1列和第2列中声明)和第二个包含多个范围的列表(在第3列,第4列和第5列中声明).我想提取第一个文件中位置落在秒文件中声明的范围内的所有行.由于位置仅在某个染色体(chr)内是唯一的,因此如果chr是相同的(即文件1中的col1与file2中的col3匹配),则必须进行测试.

档案1

chromosome position another....hundred.....columns
chr1       816 .....
chr1       991 .....
chr2       816 .....
chr2       880 .....
chr2       18768 .....
...
chr22      9736286 .....

档案2

name    identifier chromosome   start    end
GENE1   ucsc.86    chr1         800      900
GENE2   ucsc.45    chr2         700      1700
GENE3   ucsc.46    chr2         18000    19000

预期产出

chromosome position another....hundred.....columns
chr1       816 .....
chr2       816 .....
chr2       880 .....
chr2       18768 .....

我想打算做的总结(半编码):

(if $1(in file 1) matches $3(in file 2){            ##test if in the correct chr
   if ($2(in file 1) …

我有一个文件格式有点奇怪,如下所示:

Cluster 1       Score:3.96  
Category        Term        Count
GOTERM_BP_FAT   GO:0006412  34
KEGG_PATHWAY    hsa00970    9
GOTERM_BP_FAT   GO:0043038  9
GOTERM_BP_FAT   GO:0043039  9

Cluster 2       Score:3.94  
Category        Term        Count
GOTERM_BP_FAT   GO:0006414  21
KEGG_PATHWAY    hsa03010    20
GOTERM_BP_FAT   GO:0034660  16
GOTERM_BP_FAT   GO:0006399  11
GOTERM_BP_FAT   GO:0042254  10
GOTERM_BP_FAT   GO:0022613  12

...以及行之后的Cluster X行的几个"子数据帧"(包括其间的空间)和附加(此处省略)列.

我想要做的是以某种方式读取每个单独的集群,将其作为数据帧(即具有名称Category, Term, Count的数据帧),稍微操纵数据帧(主要基于计算添加列)然后写入操纵数据框架和Cluster X行到新文件的格式与它开始时的格式相同.

我已经绞尽脑汁想要做一些明智的方法,但除了分别阅读每一行并根据行的类型做不同的事情之外,我还没有真正想出任何东西,如下所示:

con  <- file('test.txt', open="r")

# Read file line for line
while ( length(currentLine <- readLines(con, n=1, warn=FALSE)) > 0 ) {
  line = strsplit(currentLine, …

我有一个基因组数据集包含基本消息,如下所示:

位置samp1 samp2 samp2 samp3 samp4 samp5 samp6 ...
posA TTTTTTT ...
posB GAAGGAA ...
posC GGGGGGG
......

该文件有100000多行,每行包含200个样本的200个碱基.
现在我想删除每个样本中具有高相似基数的位置,下面的图片是100%相同的,我将删除其中一个

我们定义了类似的比率(相似的基数)/(序列长度):

posH
CCCCCCCC posi ACCCACCC

posH和posI的相似度是6/8 = 75%根据需要,高于99%的类似比率被认为是高度相似的,并且去除了一个相似的位置.

我怎样才能有效地在python中完成这项工作？谢谢.

我使用该data.table包编写了一个脚本来解析GENCODE gtf文件的最后一列.对于那些不知道的人来说,该列包含一些键值项,每行由一个分号分隔.我正在使用的特定文件包含大约250万行.我索引前100个,然后前1000行只是为了测试脚本,输出正是我需要的.但是,尽管使用了该set功能,但运行时间并不像我预期的那么快.它是前100行的瞬间,但前1000行需要大约一两分钟.这是脚本.

#LOAD DATA.TABLE LIBRARY
require(data.table)
#READ GTF ANNOTATION FILE
info <- fread("gencodeAnnotation.gtf")
colnames(info)[9] <- "AdditionalInfo"
info <- info[1:1000]
#CREATE LIST OF 'KEYS' TO PARSE OUT
pars <- as.character(list("gene_id", "gene_type", "gene_status", "gene_name", " level ", "transcript_name", "transcript_id", "transcript_type", "transcript_support_level", "havana_gene"))
#NESTED FOR LOOP TO PARSE KEY-VALUE PAIR
for (i in 1:length(pars)) {
    for (j in 1:nrow(info)) {
        infoRow <- info[,tstrsplit(AdditionalInfo, ';', fixed = T)][j]
        headerCheck <- like(infoRow, pars[i])
        if (any(headerCheck) == TRUE) {
          keyVal <- length(tstrsplit(infoRow[[which(headerCheck …

我的python cmd脚本有问题.我不知道为什么它不起作用.我的代码可能有问题.我试图通过我的python脚本在cmdline中运行该程序.

我在bash中得到错误"sh:1:语法错误:重定向意外"

请帮助我只是生物学家:)

我使用spyder(anaconda)/ Ubuntu

#!/usr/bin/python

import sys
import os

input_ = sys.argv[1]
output_file = open(sys.argv[2],'a+')    
names = input_.rsplit('.')

for name in names:
    os.system("esearch -db pubmed -query %s  | efetch -format xml | xtract -pattern PubmedArticle  -element AbstractText >> %s" % (name, output_file))
    print("------------------------------------------")

我找到了这个perl脚本,但我有太多的序列需要分析.我想知道是否可以对其进行优化？我在它上面推出了NYTProf,看到部分"计算匹配分数","计算差距分数"和"选择最佳分数"需要花费大量时间.我是否必须修改数据结构？谢谢您的帮助.

perl脚本的参考:

    # Smith-Waterman  Algorithm
# from this website http://etutorials.org/Misc/blast/Part+II+Theory/Chapter+3.+Sequence+Alignment/3.2+Local+Alignment+Smith-Waterman/

# Smith-Waterman  Algorithm

# usage statement
die "usage: $0 <sequence 1> <sequence 2>\n" unless @ARGV == 2;

# get sequences from command line
my ($seq1, $seq2) = @ARGV;

# scoring scheme
my $MATCH    =  1; # +1 for letters that match
my $MISMATCH = -1; # -1 for letters that mismatch
my $GAP      = -1; # -1 for any gap

# initialization
my @matrix;
$matrix[0][0]{score}   = 0;
$matrix[0][0]{pointer} = "none"; …

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

例:

>chr9:128683-128744
GGATTTCTTCTTAGTTTGGATCCATTGCTGGTGAGCTAGTGGGATTTTTTGGGGGGTGTTA
>chr16:134222-134283
AGCTGGAAGCAGCGTGGGAATCACAGAATGGCCGGGACCTTAAAGGCTTTGCTTGGCCTGG
>chr16:134226-134287
GGAAGCAGCGTGGGAATCACAGAATGGACGGCCGATTAAAGGCTTTGCTTGGCCTGGATTT
>chr1:134723-134784
AAGTGATTCACCCTGCCTTTCCGACCTTCCCCAGAACAGAACACGTTGATCGTGGGCGATA
>chr16:135770-135831
GCCTGAGCAAAGGGCCTGCCCAGACAAGATTTTTTAATTGTTTAAAAACCGAATAAATGTT

此文件分为不同的部分,每个部分有2行.第一行以>(此行称为ID)开头,第二行是字母序列.我想GGAC在字母序列中搜索一个短主题(),如果它们包含主题,我想获得该行的ID(在该行之上一行).

对于上面的例子,这是预期的输出:

预期产量:

>chr16:134222-134283
>chr16:134226-134287

我试图使用以下命令在python中执行此操作,但它不返回我想要的.

infile = open('infile.txt', 'r')
ss = 'CCGA'
new = []
for line in range(len(infile)):
    if not infile[line].startswith('>'):
        for match in pattern.finder(ss):
            new.append(infile[line-1])

你知道怎么解决吗？