我是python的初學者,我不知道如何解決這個問題:
我有 60 個核苷酸的 DNA 序列以 csv 格式保存,其中每個核苷酸被 3 個二進制指示變數替換,如下所示:
A -> 1 0 0
C -> 0 1 0
G -> 0 0 1
T -> 0 0 0
所以我有一個 60*3 列和 ~3000 行的矩陣。
我想回去拿原始DNA序列。我怎樣才能做到這一點?
我嘗試了以下代碼,但它崩潰了:
mm_output = pd.DataFrame(columns=['colname_' str(i) for i in range(60)], index=range(X.shape[0]))
for i in range(X.shape[0]):
for j in np.arange(0, 180, 3):
X_sub=X.iloc[i,j:j 3]
j_tmp=int(j/3)
if np.array_equal(X_sub.values,A):
mm_output.iloc[i,j_tmp]='A'
elif np.array_equal(X_sub.values,C):
mm_output.iloc[i,j_tmp]='C'
elif np.array_equal(X_sub.values,G):
mm_output.iloc[i,j_tmp]='G'
elif np.array_equal(X_sub.values,T):
mm_output.iloc[i,j_tmp]='T'
資料集在此處收集:https : //www.openml.org/d/40670
先感謝您
uj5u.com熱心網友回復:
看起來您在 中擁有資料X,并且基于您鏈接到的檔案,X其形狀為 (3186, 180)。這是一種可以將其轉換為形狀為 (3186, 60) 的陣列的方法,其中包含字符“A”、“C”、“G”和“T”。
首先,請注意,您可以通過將每個三元組與 [1, 2, 3] 打點(即矩陣相乘)將其轉換為整數索引。例如,
In [42]: a = np.array([0, 1, 0])
In [43]: a @ [1, 2, 3] # Convert `a` to the value 2.
Out[43]: 2
In [44]: a = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]])
In [45]: a
Out[45]:
array([[1, 0, 0],
[0, 1, 0],
[0, 0, 1],
[0, 0, 0]])
In [46]: a @ [1, 2, 3] # Convert each row of `a` to a single integer.
Out[46]: array([1, 2, 3, 0])
該矩陣乘積產生的整數可以用作陣列的索引,以將它們轉換為核苷酸字符。例如
In [47]: nucs = np.array(['T', 'A', 'C', 'G'])
In [48]: nucs[[1, 2, 3, 0]]
Out[48]: array(['A', 'C', 'G', 'T'], dtype='<U1')
要將這些技術應用于您的陣列X,必須將陣列從 (3186, 180) 重新整形為 (3186, 60, 3)。(下面使用 -1reshape表示根據其他值計算出該長度。您也可以將 60 放在那里。)
In [49]: X3 = X.reshape(X.shape[0], -1, 3)
In [50]: X3.shape
Out[50]: (3186, 60, 3)
現在我們可以將矩陣乘法和索引應用于X3。創建索引陣列:
In [51]: index = X3 @ [1, 2, 3]
In [52]: index.shape
Out[52]: (3186, 60)
In [53]: index
Out[53]:
array([[2, 0, 1, ..., 1, 3, 1],
[3, 3, 0, ..., 0, 1, 2],
[3, 3, 2, ..., 2, 3, 3],
...,
[3, 3, 2, ..., 3, 2, 1],
[0, 0, 2, ..., 0, 3, 2],
[1, 3, 1, ..., 0, 3, 3]])
使用index和nucs將索引轉換為核苷酸字符:
In [54]: nucleotides = nucs[index]
In [55]: nucleotides.shape
Out[55]: (3186, 60)
In [56]: nucleotides
Out[56]:
array([['C', 'T', 'A', ..., 'A', 'G', 'A'],
['G', 'G', 'T', ..., 'T', 'A', 'C'],
['G', 'G', 'C', ..., 'C', 'G', 'G'],
...,
['G', 'G', 'C', ..., 'G', 'C', 'A'],
['T', 'T', 'C', ..., 'T', 'G', 'C'],
['A', 'G', 'A', ..., 'T', 'G', 'G']], dtype='<U1')
如果要將每行中的符號組合成單個 60 個字符的字串,可以使用
In [57]: sequences = nucleotides.view('U60')[:, 0]
In [58]: sequences.shape
Out[58]: (3186,)
In [59]: sequences
Out[59]:
array(['CTAGGCTCCAGATAGCCATAGAAGAACCAAACACTTTCTGCGTGTGTGAGAATAATCAGA',
'GGTGTTGCTCTTAGGATGTATCCCCTCAAACCTACCTGGTGGTTCTGTGCCTTCCCCTAC',
'GGCTGGACCGACCACAGCGCGTGCAGTAAGTCGGCCCCCTGCCCCGTCCTGCCCTGCCGG',
...,
'GGCCCAGCTGAGGTCTCTCTCTTCTCGCAGTTTCCATGTGGTACACACTCCCCCGCTGCA',
'TTCAGGGGAAAAAAATATCTTCTTGGCAAGGTAACACACTCTGTAAATGCATGTTCATGC',
'AGACACAGAAGTCCATCTATTACATCACTGGTGCGTTGACTCTGATTGAAGCCTTTTTGG'],
dtype='<U60')
這是沒有所有中間資料檢查的代碼:
nucs = np.array(['T', 'A', 'C', 'G'])
X3 = X.reshape(X.shape[0], -1, 3)
index = X3 @ [1, 2, 3]
nucleotides = nucs[index]
sequences = nucleotides.view('U60')[:, 0]
轉載請註明出處,本文鏈接:https://www.uj5u.com/qukuanlian/390547.html
