問題タブ [fastq]

いくつかの生物学的配列データを分析するためのワークフローを作成するために、GNUmakefile を作成しています。データは fastq と呼ばれる形式で提供され、その後、多数のクリーニングおよび分析ツールが適用されます。私が現在書いているものを添付しました。これは、クリーニング前の品質管理から、その後の品質管理までを網羅しています。私の問題は、「fastqc」コマンドを実行する方法がわからないことです。そのターゲットは、ワークフローの他のステップのいずれにも依存しないためです。

command.py、2 つのペアのファイルを一緒にマージしCA01_S1_R1.fastqますCA01_S1_R2.fastq。次に、結果を新しいディレクトリに出力し、paired.out結果のファイルに名前を付けますpaired.fastq。完全なコマンドは次のようになります

command.py -f CA01_S1_R1.fastq -r CA01_S1_R2.fastq -o paired.out

ただし、このコマンドを多くのファイルに対して実行し、すべての出力を同じディレクトリに保存したいと考えています。さらに、出力には一意の名前が必要です。したがって、これらのコマンドも効果的に実行して、ファイル 2 と 3 も送信したいと思います。

command.py -f CA02_S2_R1.fastq -r CA02_S2_R2.fastq -o paired.out

command.py -f CA03_S3_R1.fastq -r CA03_S3_R2.fastq -o paired.out

ただし、すべてのサンプルでこのコマンドをループするコードがあったとしても、すべての出力がフォルダー内paired.outにファイル名で保存されるため、コマンドは最後のペアリングの出力を上書きし続けますpaired.fastq。コマンドを介して各ファイルペアを送信し、フォルダーに入り、ファイル出力の名前を からpaired.fastqに変更しCA01_paired.fastq、すべてのファイルに対して繰り返す簡単なループを作成できますか?

次のコマンドを使用して、複数のファイルを送信できることを知っています。

このループに 2 番目の命令を追加して、このフォルダーに cd し、ファイルの名前を変更して、毎回 1 ずつインクリメントしたいと思います。インクリメント変数の設定方法がわかりません。次のようになると思います。

サブディレクトリが 94 個あり、それぞれに 1 つまたは 2 つのファイルが含まれています*.fastq。これらの各ファイルに同じ python コマンドを適用し、新しいファイルを作成する必要がありますqc_*.fastq。

各ファイルに個別にbashスクリプトを適用する方法は知っていますが、一度にすべてのファイルにコマンドを適用するbashスクリプトを書く方法があるかどうか疑問に思っています

次のタスクがあります: 30 文字の長さのパターン シーケンス (実際には DNA シーケンスであり、P30 と呼ばないようにします) から始めて、正確な P30 で始まる (^agacatacag...)、次に 29 で始まるすべての行をテキスト ファイルで検索する必要があります。 30、28、および最大 10 文字の最後の文字。パターンの最初の文字を削除して、検索を続けるだけです。簡単にするために、現在は正確な一致が必要ですが、1 つの不一致をより長く (20 ～ 30 文字の長さのパターン) 許可する方がよいでしょう。

私の現在の、かなり遅い解決策は、1 行に 1 つの切り捨てられたパターンを持つシェル ファイルを作成し、それを grep[1] することです。つまり、数 GB の巨大なテキスト ファイルを 20 倍に読み込んでおり、これには 1 日以上かかる場合があります。

Python に切り替えて、必要なすべてのパターンを含むリスト/タプルを作成し、ファイルを 1 回だけ読み取り、代わりにシーケンスごとに 20 回ループして、pypy を使用して処理を高速化できます。

• 質問 1: そのようなループよりも高速な正規表現はありますか?
• 質問 2: より高速なコンパイル済み言語に切り替えて高速化することは理にかなっていますか? （私はDlangを理解しようとしています）

[1] DNA 配列であり、検索する入力が FASTQ 形式であるため、fqgrep を使用しています: https://github.com/indraniel/fqgrep with tre ライブラリ: https://github.com/laurikari/tre/

edit_1 変更 (短縮パターン) の例。最初の数ステップ/短いパターンを示します:

または、DNA として使用する場合:

edit_2 単純な grep では実際にはうまくいきません。行番号 2 のみが一致する FASTQ 形式の各 4 行を後処理する必要があります。fqgrep を使用しない場合は、次のことを行う必要があります。

入力の 4 行を読み取ります
- 行 #2 (シーケンス) が 20 個のパターン (P30-P10) のいずれかで始まるかどうかを確認します -
一致した場合は、行 #2 と # の最初の N 文字を切り取る必要があります4、ここで、N は一致するパターンの長さを表します - 出力時に出力/ファイルへの書き込み 行 #1-$4 が一致しない場合は何もしません

社内ソリューションの場合、入力ファイルをたとえば 4M の嘘チャンクに分割し、そのように高速化する GNU 並列使用を試すことができます。しかし、新しいソフトウェアごとに他のユーザーが使用できるようにしたい場合は、エンドユーザーに広告をインストールするように依頼しています.

** edit 3 ** Vyctor の 3 つの正規表現と一致する行の簡単な例:

左から文字/DNA 塩基を削除します (または DNA の 5 プライム エンド)。これは、これらの配列が実際の酵素によって分解される方法であるためです。正規表現シーケンス自体は、いったん見つかると興味深いものではありません。望ましい出力は、正規表現の後の読み取りシーケンスです。上記の例では UPERCASE になっていますが、次のステップでゲノムにマッピングできます。このおもちゃの例とは別に、正規表現パターンの後に先験的に未知のさまざまなシーケンスが長くなっていることを強調する必要があります。現実の世界では、DNA の大文字/小文字 (すべて大文字) を処理する必要はありませんが、パターンを検索している配列の中に Ns (=未知の DNA 塩基) に遭遇する可能性があります。これらは最初の概算では無視できますが、より感度の高いバージョンのアルゴリズムでは、おそらく単純な不一致として対処する必要があります。理想的なシナリオでは、特定の位置で単純なミスマッチをカウントするのではなく、FASTQ 形式で保存された各 4 行の長いシーケンス レコードの 4 行目に保存されている DNA シーケンスの品質値を考慮して、より複雑なペナルティを計算します。http://en.wikipedia.org/wiki/FASTQ_format#Quality

しかし、これははるかに複雑であり、これまでのところ、「正規表現に完全に一致する読み取りのみを取得する」という方法で十分であり、その後の手順を分析しやすくしました。