Python實現(xiàn)將DNA序列存儲為tfr文件并讀取流程介紹
最近導(dǎo)師讓我跑模型,生物信息方向的,我一個學(xué)計算機的,好多東西都看不明白。現(xiàn)在的方向大致是,用深度學(xué)習(xí)的模型預(yù)測病毒感染人類的風(fēng)險。
既然是病毒,就需要拿到它的DNA,也就是堿基序列,然后把這些ACGT序列丟進模型里面,然后就是預(yù)測能不能感染人類,說實話,估計結(jié)果不會好,現(xiàn)在啥都是transformer,而且我看的這篇論文,我認(rèn)為僅僅從DNA序列大概預(yù)測不出什么東西。
但是就那樣吧,現(xiàn)在數(shù)據(jù)去哪里下載,需要下載什么樣的數(shù)據(jù),下載完成后怎么處理我還是一臉懵逼,但是假設(shè)上面都處理好了,然后即使把數(shù)據(jù)丟給模型,跑就完了。
也不是沒進度,目前了解到的是,我應(yīng)該使用一種叫fasta格式的文件,然后把里面的一大串ACGT序列拿出來,轉(zhuǎn)為模型可以處理的數(shù)據(jù)。然后,以后再說。
現(xiàn)在假設(shè)我已經(jīng)有了ACGT的序列,然后把它轉(zhuǎn)為模型可以處理的矩陣。
這里,我隨機生成長度為131072的基因序列,為什么是這個數(shù)字呢,因為這是之前看的 論文里的值,,暫時按照這個來做。
實現(xiàn):
首先是導(dǎo)入庫
import numpy as np import random import tensorflow as tf import inspect from typing import Any, Callable, Dict, Optional, Text, Union, Iterable import os
然后,定義一個生成長度為131072bp的函數(shù):
#隨機生成131072的dna序列
length = 131072
def randomSeq(length):
return ''.join([random.choice('ACGT') for i in range(length)])
這個函數(shù)的返回結(jié)果是長度為length的字符串,類似ACGTTGC這樣。
然后這種序列模型是沒辦法處理的,所以需要把它變成矩陣,也就用one-hot編碼。
比如ACGT這個序列,編碼成:
[ [1,0,0,0],
[0,1,0,0],
[0,0,1,0],
[0,0,0,1] ]
這樣的一個矩陣,這個就不細(xì)說了,網(wǎng)上很多資料。
然后,我從別人的代碼中抄了一個函數(shù),很好用。
#DNA序列轉(zhuǎn)為one-hot編碼,可以直接拿來用
def one_hot_encode(sequence: str,
alphabet: str = 'ACGT',
neutral_alphabet: str = 'N',
neutral_value: Any = 0,
dtype=np.float32) -> np.ndarray:
"""One-hot encode sequence."""
def to_uint8(string):
return np.frombuffer(string.encode('ascii'), dtype=np.uint8)
hash_table = np.zeros((np.iinfo(np.uint8).max, len(alphabet)), dtype=dtype)
hash_table[to_uint8(alphabet)] = np.eye(len(alphabet), dtype=dtype)
hash_table[to_uint8(neutral_alphabet)] = neutral_value
hash_table = hash_table.astype(dtype)
return hash_table[to_uint8(sequence)]
這是一個嵌套函數(shù)了,仔細(xì)研究下還是可以理解的,我就不說了,會用就行了。
簡單講一下參數(shù)的意思:
sequence:字符串類型,就是輸入的堿基序列。
alphabet: str = ‘ACGT’ :詞表,一共只需要這四個詞
neutral_alphabet: str = ‘N’,
neutral_value: Any = 0,
上面這兩一起用,就是說遇到N這個堿基就會編碼成[0,0,0,0]的向量。
dtype=np.float32,這個就是內(nèi)部元素值的類型。
簡單生成一下:

然后輸入序列長度是131072bp,所以輸入的矩陣就是131072x4的矩陣,現(xiàn)在來把序列變?yōu)榫仃嚒?/p>
編碼成one-hot矩陣
dnaVec = one_hot_encode(dna)

現(xiàn)在DNA序列已經(jīng)變成了矩陣,接下來需要把這一條序列,也就是一個樣本數(shù)據(jù),變成TensorFlow中的TFRecord文件格式。TFRecord 是 TensorFlow 中的數(shù)據(jù)集存儲格式。當(dāng)我們將數(shù)據(jù)集整理成 TFRecord 格式后,TensorFlow 就可以高效地讀取和處理這些數(shù)據(jù)集,從而幫助我們更高效地進行大規(guī)模的模型訓(xùn)練。
關(guān)于tfr文件的處理,我就不在細(xì)說了,總之現(xiàn)在我們需要構(gòu)建example。
在此之前,我們需要先這么做:
#給出結(jié)果的tfr文件的路徑 path = '/content/drive/MyDrive/test_Enformer/result.tfr' #dna的numpy數(shù)組轉(zhuǎn)成字節(jié)流,這樣才能存儲 dnaVec = dnaVec.tobytes()
接下來就是把這個字節(jié)流數(shù)據(jù)寫入到tfr文件中,這里同時寫入這條數(shù)據(jù)的label中,我的問題是給一個Dna序列,預(yù)測是或者不是的二分類問題,所以我同時把這條dna序列對應(yīng)的真實標(biāo)簽也寫進去,但是我是隨機從0,1中選擇一個。
from tensorflow.core.example.feature_pb2 import BytesList
with tf.io.TFRecordWriter(path) as writer:
feature = {
#序列使用的是tf.train.BytesList類型
'sequence':tf.train.Feature(bytes_list=tf.train.BytesList(value=[dnaVec])),
#label是隨機生成的0,或者1
'label':tf.train.Feature(int64_list=tf.train.Int64List(value=[np.random.choice([0,1])]))
}
example = tf.train.Example(features=tf.train.Features(feature=feature))
writer.write(example.SerializeToString())這部分的代碼執(zhí)行結(jié)束后,就已經(jīng)把dna序列以及對應(yīng)的標(biāo)簽寫入了tfr文件中,不過這個tfr文件中只有一個example,你可以寫更多個。
剛剛寫入的tfr文件

到這里,相當(dāng)于已經(jīng)把數(shù)據(jù)準(zhǔn)備好了,接下來就是讀取數(shù)據(jù)。
#從剛才的路徑中加載數(shù)據(jù)集 dataset = tf.data.TFRecordDataset(path)
#定義Feature結(jié)構(gòu),告訴解碼器每個Feature的類型是什么
feature_description = {"sequence": tf.io.FixedLenFeature((), tf.string),
"label": tf.io.FixedLenFeature((), tf.int64)}
#將 TFRecord 文件中的每一個序列化的 tf.train.Example 解碼
def parse_example(example_string):
#解析之后得到的example
example = tf.io.parse_single_example(example_string,feature_description)
#example['sequence']還是字節(jié)流的形式,重新轉(zhuǎn)為數(shù)字向量
sequence = tf.io.decode_raw(example['sequence'], tf.float32)
sequence = tf.reshape(sequence,(length,4)) #形狀需要重塑,不然就是一個長向量
label = tf.cast(example['label'],tf.int64) #標(biāo)簽對應(yīng)的類型轉(zhuǎn)換
#每一天example解析后返回對應(yīng)的一個字典
return {
'sequence':sequence,
'label': label
}#把parse_example函數(shù)映射到dataset中的每個example, #這里的dataset中只有一個example dataset = dataset.map(parse_example)
此時的dataset是一個可以遍歷的對象,內(nèi)部元素可以認(rèn)為是解析完成后的example。

這個字典有兩個鍵sequence和lable,對應(yīng)著序列矩陣和標(biāo)簽值

這就是可以用來訓(xùn)練的數(shù)據(jù)。
到此這篇關(guān)于Python實現(xiàn)將DNA序列存儲為tfr文件并讀取流程介紹的文章就介紹到這了,更多相關(guān)Python存儲tfr文件內(nèi)容請搜索腳本之家以前的文章或繼續(xù)瀏覽下面的相關(guān)文章希望大家以后多多支持腳本之家!
相關(guān)文章
python爬蟲入門教程--快速理解HTTP協(xié)議(一)
http協(xié)議是互聯(lián)網(wǎng)里面最重要,最基礎(chǔ)的協(xié)議之一,我們的爬蟲需要經(jīng)常和http協(xié)議打交道。下面這篇文章主要給大家介紹了關(guān)于python爬蟲入門之快速理解HTTP協(xié)議的相關(guān)資料,文中介紹的非常詳細(xì),需要的朋友可以參考借鑒,下面來一起看看吧。2017-05-05
詳解如何基于Pyecharts繪制常見的直角坐標(biāo)系圖表
pyecharts是基于前端可視化框架echarts的Python可視化庫,下面這篇文章主要給大家介紹了關(guān)于如何基于Pyecharts繪制常見的直角坐標(biāo)系圖表的相關(guān)資料,文中通過實例代碼介紹的非常詳細(xì),需要的朋友可以參考下2022-04-04
在Python運行時動態(tài)查看進程內(nèi)部信息的方法
今天小編就為大家分享一篇在Python運行時動態(tài)查看進程內(nèi)部信息的方法,具有很好的參考價值,希望對大家有所幫助。一起跟隨小編過來看看吧2019-02-02
Python語法學(xué)習(xí)之線程的創(chuàng)建與常用方法詳解
本文主要介紹了線程的使用,線程是利用進程的資源來執(zhí)行業(yè)務(wù),并且通過創(chuàng)建多個線程,對于資源的消耗相對來說會比較低,今天就來看一看線程的使用方法具體有哪些吧2022-04-04
Python基于ThreadingTCPServer創(chuàng)建多線程代理的方法示例
這篇文章主要介紹了Python基于ThreadingTCPServer創(chuàng)建多線程代理的方法,結(jié)合實例形式分析了Python使用ThreadingTCPServer模塊實現(xiàn)多線程代理功能進行網(wǎng)絡(luò)請求響應(yīng)的相關(guān)操作技巧,需要的朋友可以參考下2018-01-01
跟老齊學(xué)Python之永遠(yuǎn)強大的函數(shù)
Python程序中的語句都會組織成函數(shù)的形式。通俗地說,函數(shù)就是完成特定功能的一個語句組,這組語句可以作為一個單位使用,并且給它取一個名字,這樣,我們就可以通過函數(shù)名在程序的不同地方多次執(zhí)行(這通常叫做函數(shù)調(diào)用),卻不需要在所有地方都重復(fù)編寫這些語句。2014-09-09

