-3

私は2つのファイルを持っています.1つはfastaファイルで、もう1つはfastqファイルです。fastaを取得し、各行を読み取り、fastqファイル内の各行を検索して、一番上の行と一番下の行を出力したいと考えています。これは私が持っているものです

fasta ファイル

読み取り1

あああああああああああああああああああああああああああああああああああああ

ああああああああああああああああああああああああああああああああああああああ

ああああああああああああああああああああああああああああああああああああああ

あああああああああああああああああああああああああああああああああああああ

あああああああああああああああああああああああああああああああああああああああ

for seq in `cat sequences`;do grep -A 2 -B 1 $seq FAP.1.txt;done

@DH1DQQN1:269:C1UKCACXX:1:1107:20386:6577 1:N:0:TTAGGC

ああああああああああああああああああああああああああああああああああああああ

+

CCCFFFFFHGHHHJIJHFDDDB173@8815BDDB################

@DH1DQQN1:269:C1UKCACXX:1:1114:5718:53821 1:N:0:TTAGGC

あああああああああああああああああああああああああああああああああああああ

+ ;@?DBD<@@FFHHH<

@DH1DQQN1:269:C1UKCACXX:1:1209:10703:35361 1:N:0:TTAGGC

あああああああああああああああああああああああああああああああああああああ

+

@@@FFFFFHGHHHGIJHFDDDDDBDD69@6B-707537BDDDB75@@85

@DH1DQQN1:269:C1UKCACXX:1:1210:18926:75163 1:N:0:TTAGGC

ああああああああああああああああああああああああああああああああああああああ

@CFFFFFFHHHHHJJJHFDDD@77BDDDDB077007@B############

このことから、それが 2 回表示されていることがわかりますがAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA、印刷したいのは 1 回だけです。どうやってやるの?

最終出力ファイル

@DH1DQQN1:269:C1UKCACXX:1:1107:20386:6577 1:N:0:TTAGGC

ああああああああああああああああああああああああああああああああああああああ

+

CCCFFFFFHGHHHJIJHFDDDB173@8815BDDB################

@DH1DQQN1:269:C1UKCACXX:1:1114:5718:53821 1:N:0:TTAGGC

あああああああああああああああああああああああああああああああああああああ

+

;@?DBD<@@FFHHH<

@DH1DQQN1:269:C1UKCACXX:1:1210:18926:75163 1:N:0:TTAGGC

ああああああああああああああああああああああああああああああああああああああ

+

@CFFFFFFHHHHHJJJHFDD@77BDDDDB077007@B

4

3 に答える 3

0

これを自分で行う強い理由がない限り、Biopythonを使用してください。

ファスト:

AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAC
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAG
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAGA
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAGG

fastq (あなたのものに基づいていますが、出力の形式が悪いため同一ではありません):

@DH1DQQN1:269:C1UKCACXX:1:1107:20386:6577 1:N:0:TTAGGC
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAC
+
CCCFFFFFHGHHHJIJHFDDDB173@8815BDDB###############
@DH1DQQN1:269:C1UKCACXX:1:1114:5718:53821 1:N:0:TTAGGC
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
+
CCCFFFFFHGHHHJIJHFDDDB173@8815BDDB###############
@DH1DQQN1:269:C1UKCACXX:1:1209:10703:35361 1:N:0:TTAGGC
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
+
@@@FFFFFHGHHHGIJHFDDDDDBDD69@6B-707537BDDDB75@@85
@DH1DQQN1:269:C1UKCACXX:1:1210:18926:75163 1:N:0:TTAGGC
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAG
+
@CCFFFFFHHHHHJJJHFDDD@77BDDDDB077007@B###########

コード:

from Bio import SeqIO

with open("fasta") as fh:
    fasta = fh.read().splitlines()

seen = set()

for record in SeqIO.parse(open('fastq'), 'fastq'):
    seq = str(record.seq)
    if seq in fasta and seq not in seen:
        seen.add(seq)
        print record.format('fastq')

EDIT : 上記は、fasta ファイルではなく、fastq ファイルの順序でレコードを出力します。順序が重要でない場合は、その方法を使用する必要があります。それ以外の場合は、キーが FASTA ファイル内のインデックスであるディクショナリにレコードを追加し、最後にそれらすべてを出力して、ディクショナリをソートできます。

from Bio import SeqIO
import sys

with open("fasta") as fh:
    fasta = fh.read().splitlines()

seen = set()
records = {}

for record in SeqIO.parse(open('fastq'), 'fastq'):
    seq = str(record.seq)
    if seq in fasta and seq not in seen:
        seen.add(seq)
        records[fasta.index(seq)] = record

for record in sorted(records):
    sys.stdout.write(records[record].format('fastq'))

(ここでは、余分な改行を避けるために、のsys.stdout.write代わりに, も使用しています。)print

于 2013-08-02T23:14:58.603 に答える
0

一意のシーケンスのみを含む fastq が必要なようですね。

これは恐ろしく非効率的な方法ですが、うまくいくはずです。fastq ファイルをリストとして保存するので、大きすぎないことを願っています。品質スコアなどではなく、重複シーケンスのみをスローします。

fastqFile = list(open(fastq))
out = []
output = open('output.fastq', 'at')

for lineNum, line in enumerate(fastqFile):
    if lineNum < 4:
        out.append(line)
        output.write(line)
    else:
        if line not in out and lineNum % 4 != 3:
            output.write(fastqFile[lineNum - 1])
            output.write(line)
            output.write(fastqFile[lineNum + 1])
            output.write(fastqFile[lineNum + 2])
            out.append(fastqFile[lineNum - 1])
            out.append(line)
            out.append(fastqFile[lineNum + 1])
            out.append(fastqFile[lineNum + 2])
于 2013-08-02T21:59:59.743 に答える