【问题标题】:How to efficiently iterate through a large string and record the index of every 3-letter substring?如何有效地遍历一个大字符串并记录每个 3 字母子字符串的索引?
【发布时间】:2020-04-16 21:05:37
【问题描述】:

我有一个非常大的序列,我以字符串的形式读入(250,000,000 个字母)。字母是 G、A、C、T。

例如:

'GACTCGTAGCTAGCTG'

我想创建一些方法来存储序列中每个 3 个字母子字符串的索引,以便稍后在另一个函数中使用。

例如:

{'GAC': 1, 'ACT': 2 'CTC':3, 'TCG':4 ...}

对于我目前的方法,我的问题是我还没有找到一种有效的方法来存储序列中每个 3 个字母子字符串的索引。一旦我知道每个子字符串的索引,我将根据我拥有的给定概率随机选择其中一些并将它们更改为另一个已知的子字符串。

我尝试使用 for 循环进行迭代,在滑动窗口中一次使用一个 3 个字母的子字符串,并将位置保存为字典值,并以 3 个字母的子字符串作为键。另外,当我保存这个字典文件时,我一直在使用 pd.to_csv 但它似乎效率很低,并创建了一个 8 GB 的文件。然而,对于一个非常大的字符串(250,000,000 个字母)来说,这需要很长时间。

【问题讨论】:

  • 我看起来你正在制作一个字典,其中键是三个字母的子字符串,但值是简单的数字。当子字符串出现在您的示例中的多个位置(例如AGC)时,您会怎么做。
  • 你可以在列表中做相反的事情... ['GAC', 'ACT', 'CTC': 'TCG':...] 所以索引 0 是 gac,1 是 act ...所以它是相同的但相反...如果你必须保持你的格式,你需要一个列表作为每个索引的值
  • stackoverflow.com/questions/39944594/… 也许问题不在于python,而是数据太大而无法开始?我还建议您使用子字符串的位置作为键而不是值。

标签: python string


【解决方案1】:

您可以将长字符串拆分为 3 个字母的子字符串:

string='GACTCGTAGCTAGCT'
substrings=[string[3*x:3*x+3] for x in range(int(len(string)/3))]

substrings 将是:

['GAC', 'TCG', 'TAG', 'CTA', 'GCT']

将这些索引添加到另一个列表中:

indices=[x for x in range(int(len(string)/3))]

这只会简单地产生:

[0, 1, 2, 3, 4]

子字符串列表中的第n个元素将对应索引列表中的第n个元素。

关于如何将文件放入字符串变量,您可能想查看: How to read a text file into a string variable and strip newlines?

【讨论】:

    【解决方案2】:

    您的示例输出表明每个子序列只需要一个索引,但您的描述暗示您可能需要所有索引。

    这是一个函数,它将构建一个字典,其中包含存在的子序列(给定长度)的所有索引:

    from collections import defaultdict
    def sequenceDict(seq,count):
        result = defaultdict(list)
        for i,subseq in enumerate(zip(*(seq[i:] for i in range(count)))):
            result["".join(subseq)].append(i)
        return result
    
    r = sequenceDict('GACTCGTAGCTAGCTG',3)
    print(r)
    
    # {'GAC': [0], 'ACT': [1], 'CTC': [2], 'TCG': [3], 'CGT': [4], 'GTA': [5], 'TAG': [6, 10], 'AGC': [7, 11], 'GCT': [8, 12], 'CTA': [9], 'CTG': [13]})
    

    如果你真的只想要每个 3 字母子序列的第一个索引,那么使用单个字典理解可以更快地获得字典:

    from itertools import product
    { ss:sequence.index(ss) for p in product(*["ACGT"]*3)for ss in ["".join(p)] if ss in sequence}
    

    我对 2.5 亿个字母的随机序列进行了性能测试,可以在几微秒内获得单索引字典。获取所有索引需要一分钟多一点(使用上面的函数):

    import time
    
    size = 250_000_000
    print("loading sequence...",size)
    start = time.time()
    import random
    sequence = "".join(random.choice("ACGT") for _ in range(size))
    print("sequence ready",time.time()-start)
    
    
    start = time.time()
    from itertools import product
    seqDict = { ss:sequence.index(ss) for p in product(*["ACGT"]*3)for ss in ["".join(p)] if ss in sequence}
    print("\n1st index",time.time()-start)
    
    start = time.time()
    r = sequenceDict(sequence,3)
    print("\nall indexes",time.time()-start)
    

    输出:

    loading sequence... 250000000
    sequence ready 193.82172107696533
    
    1st index 0.000141143798828125
    
    all indexes 71.74848103523254
    

    鉴于加载序列的时间比构建索引的时间长得多,您可能会放弃存储该索引字典并每次都从源数据中重建它(您似乎仍然需要为您的进程加载)

    您也可以只存储一个计数字典并根据需要提取索引:

    此函数获取每个子序列的出现次数:

    from collections import Counter
    def countSubSeqs(seq,size):
        return Counter("".join(s) for s in zip(*(seq[i:] for i in range(size))))
    

    它的运行时间与 sequenceDict 函数大致相同,但生成的字典要小得多。

    要获取特定子序列(包括重叠位置)的索引,您可以使用:

    subSeq  = "ACT"
    indexes = [ i for i in range(len(sequence)) if sequence[i:i+3]==subSeq ] 
    

    如果您不需要立即为所有子序列创建所有索引,则可以相应地构建代码并仅在需要时获取索引(并可能将它们存储在字典中以供查询和重用)

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 2021-02-04
      • 1970-01-01
      • 2014-04-02
      • 2013-01-18
      • 2021-06-15
      • 1970-01-01
      • 2010-09-09
      相关资源
      最近更新 更多