时间戳

首页 > 解决方案 > 在矩阵中管理和子串多个字符串（Python & Numpy）

python - 在矩阵中管理和子串多个字符串（Python & Numpy）

问题描述

我有一个名为“test.fa”的文件，内容如下：

>Sequence 1
TCAGAACCAGTTATAAATTTATCATTTCCTTCTCCACTCCT
>Sequence 2
CCCACGCAGCCGCCCTCCTCCCCGGTCACTGACTGGTCCTG
>Sequence 3
TCGACCCTCTGGAACCTATCAGGGACCACAGTCAGCCAGGCAAG
>Sequence 4
AAAACACTTGAGGGAGCAGATAACTGGGCCAACCATGACTC

这个测试文件只有 8 行，但可以有更多。我现在也不知道所有序列的长度。所以首先我计算行数，并根据行数创建一个矩阵。

import numpy as np
filename = "test.fa"
with open(filename, 'r') as file:
    n1 = 0
    n2 = 0
    for line in file.readlines():
        if line[0] != '>':
            m1 = 1
            n1 = n1 + m1
        else:
            m2 = 1
            n2 = n2 + m2
    seq = np.chararray((n1,2),itemsize = 99)

接下来，我将值添加到矩阵中。

with open(filename, 'r') as file:
    n3 = -1
    n4 = -1
    for line in file.readlines():
        if line[0] != '>':
            m3 = 1
            n3 = n3 + m3
            seq[n3,1] = line
        else:
            m4 = 1
            n4 = n4 + m4
            seq[n4,0] = line

如果我调用'seq'，我会得到：

chararray([[b'>Sequence 1', b'TCAGAACCAGTTATAAATTTATCATTTCCTTCTCCACTCCT'],
       [b'>Sequence 2', b'CCCACGCAGCCGCCCTCCTCCCCGGTCACTGACTGGTCCTG'],
       [b'>Sequence 3', b'TCGACCCTCTGGAACCTATCAGGGACCACAGTCAGCCAGGCAAG'],
       [b'>Sequence 4', b'AAAACACTTGAGGGAGCAGATAACTGGGCCAACCATGACTC']], dtype='|S99')

seq[2,1]: b'TCGACCCTCTGGAACCTATCAGGGACCACAGTCAGCCAGGCAAG'
seq[2,1][0] : 84
seq[2,1][0:8] : b'TCGACCCT'
seq[2,8][8] : 67
seq[2,1][8:9] : C

这感觉不是管理这些数据的最佳方式（为什么 [2,1][0] 返回 84？）。我不确定项目大小为 99 的 np.chararray 是最好的方法。

我的问题是：组织/管理这些数据的更好方法是什么？我最终需要计算每个核苷酸的出现次数（即有多少 As、Ts、Cs、Gs），并从每个序列中提取子串。就上下文而言，这与期望最大化和吉布斯抽样有关。

标签： pythonnumpysubstring

解决方案

它返回 84，而类型是 S99，因为找到的最大长度是 99，所以类型是 S99，你从文件中读取的方式不是那么正确，管理你的数据的更好方法是使用dict，我不知道是否它将满足您的下一个需求，但这里是：

seq={}
with open("file") as f:
    for line in f:
        seq[line]=next(f)

这里它将一次读取两行，所以第一行line有序列 1，next(f) 有 seq TGC .....，如果发生文件的顺序不一致，请使用 if 语句

推荐阅读