Python入門のファイル操作(一)

前回はpythonを初歩的に認識しましたが、懸念を残しました.生物情報の中で最初に身につけなければならない基礎技能は、間違いなく、この文書の操作です.
本稿では、生物情報の基礎となるファイルの読み書き、後続の補足ファイル権限管理、ファイルの移動、複製などの操作(まず餅を描き、穴を埋める暇がある)について簡単に紹介します.
ファイルの読み書き
日常の生物情報学分析では、シーケンシングデータfastqファイル、例えば結果bamまたはsamファイル、gtfまたはgff 3フォーマットの注釈ファイル、さらにはデータベースファイルなど、いくつかのファイルに直接直面して操作します.
これらのファイルに直面するとき、インタラクティブなプログラミングは力不足で、多くの場合、スクリプトプログラミングを選択します.(インタラクティブプログラミングとスクリプトプログラミングの違いについては別の文章を参照してください.もちろんまだ書いていません!(lll￢ω￢)).
スクリプトプログラミングでまず遭遇するのはファイルの読み書きであり,次にパラメータの伝達である.
一般的なファイルはopenで開くことができ、readlineで読み取ることができます.

#    
f = open('/path/to/file','r') #       
line = f.readline() #     
linse = f.readlines() #      
f.close() #        

#    
f = open('/path/to/file','w') #      
f.write(string) #   string      
f.close() #        

ここではwith文を推奨します.コード量を減らすだけでなく、いくつかの異常を処理することもできます.

with open('/path/to/file','r') as f:
    line = f.readline() #     
    linse = f.readlines() #      

生物情報の一般的な主流のファイルフォーマットはfasta,fastq,sra,Genebank,gtf,gff 3,sam,bamなどである.
fastaファイル読み書き
fastaファイルは一般的に>で始まり、第1行は一般的に配列名、またはその他の情報などであり、第2行はその塩基配列で始まる.次のようになります.

# test.fa
>MSTR.8.1 gene=MSTR.8
CCACAATTAATCAGAGGTAGCGACTACCCTATTTCTTCATTTCTTATGATTCTTACTATTTAAGCACTCA
TTTTCCATACCCCCACCTTTTACTTTCCTCTATTTTCCTCTCATCGTTAACTGATTCTGTTCAAGTTATA
GGGTTGGTTCTCAATCCTTGATATGTTGCTGGAATCCCTTTGGATTCGCTAAAAGGAGGATTTCAAGGGT
TTTCCAGAGTTTCCGACATTTGATTAAGGCGATCGGTGGATTTGAACTTGATTTCGATTGTTTCACAGTA
CCACCATGTATGGGTAATCGAGGATTTTAGGGAGGAGAGATTTGAGCCTCAAGGAGGAAGATAACCAGTT
GGTATTTCATGGTAGTAGAATCTTCAGCAATTTCTAGGTATGTTTAGAGGGT
>MSTR.26.1 gene=MSTR.26
CACCAGATCAAGATGACAGTAACTACAATTTAAGCAGATCGAACAGGAAACATGTAGCCATCTTTCATGT
ATGAAAGGTCCAAATACTAGACTCTACTTGTGCATCCCTCAGACTACATTTATTTATCTTGCCGAAATGT
TAAAAATGAAGTTTCAGTTTCTCATGTAAAGTTTAGTTTCAGTTTCATTGAGATGAATCCATGGATGTAG
CTAGTCATCCATTCTCGGGAAACAACAACCACC

fastaファイルの読み取り方法には、一般的に次の2つの方法があります.

は、fastaファイルの各シーケンスを辞書として格納し、辞書のキーはfastaシーケンスの名前であり、値はそのアルカリ配列である.

#!/usr/bin/evn python
# -*- coding: UTF-8 -*-
with open("test.fa",'r') as fa:
    sequences = {}
    for line in fa:
        if line.startswith(">"):
            name = line.rstrip("
")
            sequences[name] = ""
        else:
            sequences[name] = sequences[name] + line.rstrip("
")
    print (sequences)

出力は次のとおりです.

{'>MSTR.26.1 gene=MSTR.26': 'CACCAGATCAAGATGACAGTAACTACAATTTAAGCAGATCGAACAGGAAACATGTAGCCATCTTTCATGTATGAAAGGTCCAAATACTAGACTCTACTTGTGCATCCCTCAGACTACATTTATTTATCTTGCCGAAATGTTAAAAATGAAGTTTCAGTTTCTCATGTAAAGTTTAGTTTCAGTTTCATTGAGATGAATCCATGGATGTAGCTAGTCATCCATTCTCGGGAAACAACAACCACC', '>MSTR.8.1 gene=MSTR.8': 'CCACAATTAATCAGAGGTAGCGACTACCCTATTTCTTCATTTCTTATGATTCTTACTATTTAAGCACTCATTTTCCATACCCCCACCTTTTACTTTCCTCTATTTTCCTCTCATCGTTAACTGATTCTGTTCAAGTTATAGGGTTGGTTCTCAATCCTTGATATGTTGCTGGAATCCCTTTGGATTCGCTAAAAGGAGGATTTCAAGGGTTTTCCAGAGTTTCCGACATTTGATTAAGGCGATCGGTGGATTTGAACTTGATTTCGATTGTTTCACAGTACCACCATGTATGGGTAATCGAGGATTTTAGGGAGGAGAGATTTGAGCCTCAAGGAGGAAGATAACCAGTTGGTATTTCATGGTAGTAGAATCTTCAGCAATTTCTAGGTATGTTTAGAGGGT'}

サードパーティ製モジュールを使用Biopython SeqIOモジュールを使用していますが、あまり使用しない前にインストールする必要があるようですが、どのようにインストールするかは自分でgoogleしてください.

#!/usr/bin/evn python
# -*- coding: UTF-8 -*-
from Bio import SeqIO

with open ("test.fa",'r') as fa:
    for record in SeqIO.parse(fa,'fasta'):
        print(record)

出力は次のとおりです.

ID: MSTR.8.1
Name: MSTR.8.1
Description: MSTR.8.1 gene=MSTR.8
Number of features: 0
Seq('CCACAATTAATCAGAGGTAGCGACTACCCTATTTCTTCATTTCTTATGATTCTT...GGT', SingleLetterAlphabet())
ID: MSTR.26.1
Name: MSTR.26.1
Description: MSTR.26.1 gene=MSTR.26
Number of features: 0
Seq('CACCAGATCAAGATGACAGTAACTACAATTTAAGCAGATCGAACAGGAAACATG...ACC', SingleLetterAlphabet())

正しく読み取れば、他の必要な操作ができます.
fastqファイル読み書き
fastqファイルフォーマット:各シーケンスは4行あり、第1行は一般的に@記号で始まり、主にシーケンシングのいくつかの情報はシーケンス名として簡単に理解することができ、第2行はその塩基配列であり、第3行も名称であり、一般的に空間を節約するために+号で表され、第4行はシーケンスの各塩基のシーケンシング品質値であり、第2行の塩基と1つに対応している.

@FCD1VGJACXX:2:1101:5762:2000#AAGTCGCG/1
NGGGAGACCGGGGTTCGATTCCCCGACGGGTCGTATGCCGTCTTCTGCT
+
BS\ceeeegggggffgighiiiihhUegiigggeeeeddcccccccccc
@FCD1VGJACXX:2:1101:7844:2000#AAGTCGCG/1
NTGGCTTCCTAAGCCAGGGATTGTGGGTTCGTATGCCGTCTTCTGCTTG
+
BPYa^accgeeggihfhhhfhiQQQQPPPPPP^XYefgeeWaeghhgih

fastqファイルについては、ほとんどの場合、主にアルカリ基の品質を考慮しています.通常、その配列についてアルカリ基統計、品質統計などを行います.fastqファイルについても、同じように2つの方法で読み取ることができます.
1.一般的な方法

#!/usr/bin/evn python
# -*- coding: UTF-8 -*-
# style 1
with open("test.fq",'r') as fq:
    lines = fq.readlines()
    head = [item[:-1] for item in lines[::4]]
    read = [item[:-1] for item in lines[1::4]]
    qual = [item[:-1] for item in lines[3::4]]
    
    sequence = dict(zip(read, qual))
    print(sequence)
    
# style 2
with open('test.fq','r') as fq:
    sequence = {}
    name = None
    for i,line in enumerate(fq):
        if i % 4 == 0:
            name = line.rstrip('
')
        if i % 4 == 1:
            read = line.rstrip('
')
        if i % 4 == 3:
            sequence[(name,read)] = line.rstrip('
')
    for key in sequence:
        print (key[0],key[1],sequence[key])

出力は次のとおりです.

# style 1
{'NTGGCTTCCTAAGCCAGGGATTGTGGGTTCGTATGCCGTCTTCTGCTTG': 'BPYa^accgeeggihfhhhfhiQQQQPPPPPP^XYefgeeWaeghhgih', 'NGGGAGACCGGGGTTCGATTCCCCGACGGGTCGTATGCCGTCTTCTGCT': 'BS\\ceeeegggggffgighiiiihhUegiigggeeeeddcccccccccc'}
# style 2
('@FCD1VGJACXX:2:1101:5762:2000#AAGTCGCG/1', 'NGGGAGACCGGGGTTCGATTCCCCGACGGGTCGTATGCCGTCTTCTGCT', 'BS\\ceeeegggggffgighiiiihhUegiigggeeeeddcccccccccc')
('@FCD1VGJACXX:2:1101:7844:2000#AAGTCGCG/1', 'NTGGCTTCCTAAGCCAGGGATTGTGGGTTCGTATGCCGTCTTCTGCTTG', 'BPYa^accgeeggihfhhhfhiQQQQPPPPPP^XYefgeeWaeghhgih')

2.サードパーティモジュール読み出し
fastaファイルを読み込むモジュールと同様にBiopythonを使用

#!/usr/bin/evn python
# -*- coding: UTF-8 -*-
from Bio import SeqIO

with open ("test.fq",'r') as fq:
    for record in SeqIO.parse(fq,'fastq'):
        print(record)

出力は次のとおりです.

ID: FCD1VGJACXX:2:1101:5762:2000#AAGTCGCG/1
Name: FCD1VGJACXX:2:1101:5762:2000#AAGTCGCG/1
Description: FCD1VGJACXX:2:1101:5762:2000#AAGTCGCG/1
Number of features: 0
Per letter annotation for: phred_quality
Seq('NGGGAGACCGGGGTTCGATTCCCCGACGGGTCGTATGCCGTCTTCTGCT', SingleLetterAlphabet())
ID: FCD1VGJACXX:2:1101:7844:2000#AAGTCGCG/1
Name: FCD1VGJACXX:2:1101:7844:2000#AAGTCGCG/1
Description: FCD1VGJACXX:2:1101:7844:2000#AAGTCGCG/1
Number of features: 0
Per letter annotation for: phred_quality
Seq('NTGGCTTCCTAAGCCAGGGATTGTGGGTTCGTATGCCGTCTTCTGCTTG', SingleLetterAlphabet())

sam/bamファイルの読み取りについては、サードパーティモジュール:pysamを使用して、具体的な練習は後で話しています.
gtfおよびgff 3注釈情報ファイルの読み出しは、一般テキストファイルの読み出し方式と一致する.