2017-12-31 260 views
0

ファイル1とファイル2の2つのファイルからsamtoolsを使用してパイルアップファイルを作成しようとしています。BASH別のファイルのある列からパイプされた値を使用してパイルアップファイルを再帰的に作成する

Iの形式次の名前付き44個のファイルを有し、その結果、染色体によりファイル1とfile2を分割しました:

chr${c}.${TISSUE}_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY 

$ {C}は1と22の間の数であり、そして$組織は結腸のいずれかでありますまたは筋肉 - 結腸の22の染色体、および筋肉の22の染色体。私は; chr1.colon_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY chr2.colon_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY

. 
. 
. 

chr22.colon_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY 
chr1.muscle_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY 
. 
. 
. 

これらのファイルは、2つの列で構成、最初はちょうど染色体数を示し、第2列ですその染色体上の位置。私は;

(例えば、 "chr2.colon_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY"用)ファイルの各行について
head chr2.colon_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY 
chr2 103977 
chr2 112051 
chr2 126199 
chr2 146288 
chr2 147797 
chr2 147822 
chr2 148548 
chr2 148525 
chr2 158189 
chr2 158188 

、Iは、位置を取るカラム2から、「X」を呼び出し、そしてa-bの範囲を得るためにそれを使用する必要があり、a=x-5及びb=x+5。例えば、

samtools mpileup -f [REFERENCE GENOME] File1 File2 -r chr${c}:a-b 

Iは第2染色体、位置103977(上の行1)で探していたとします。私は、次のスクリプトにこれらの値を接続します。その後、私のスクリプトは

です。基本的にループ内のループ内のループです。何かのように、

for t in $(colon, muscle) 
do 
    for c in $seq (1 22) 
    do 
    for item (or maybe row?) in 
     chr${c}.${t}_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY 
    do 
     awk '{print $2}' | something something something 
     x= position in col 2, a=x-5 b=x+5 
     samtools mpileup -f [REFERENCE GENOME] File1 File2 -r chr${c}:a-b 
    done 
    done 
done 
... 

ありがとうございます。私はLinuxでの作業にはまったく新しいので、基本的にコンピュータサイエンストレーニングはありません。

+2

こんにちは、編集して、ポストエディタのコード(中かっこ)機能を使用して読みやすくしてください。 質問は読めません。整理してください。 他の質問を見て、適切な質問を書く方法を学ぶことをお勧めします。 GL :) – Blacky

答えて

1

awkは一度にラインを処理し、私はつまり

for t in colon muscle; do 
    for c in $(seq 1 22); do 
     awk '{ print $2-5 "-" $2+5 }' chr${c}.${t}_run1_en2hic_PE1.bam.sorted.bam.breaks_COL1_and_COL2_ONLY | 
     while read -r range; do 
      samtools mpileup -f [REFERENCE GENOME] File1 File2 -r chr${c}:$range 
     done 
    done 
done 

のようなもののために行くだろう、Awkのは、ファイル全体を処理し、最終的while read -r rangeループに一度出力の1行を供給します。

最初にこれらのファイルをどのように分割するか分かりませんが、代わりにFile1File2に直接作業した場合、これはかなり簡単になる可能性があります。

また、おそらく外側のループを避けて、すべての*_ONLYファイルでAwkを直接実行することもできます。現在のファイル名はAwkの内部変数FILENAMEから得ることができますが、この場合は最初のフィールドを使用するだけです。あなたが直接$1を使用できない場合

awk '{ print $1 ":" $2-5 "-" $2+5 }' *_ONLY | 
while read -r chrrange; do 
    samtools mpileup -f [REFERENCE GENOME] File1 File2 -r "$chrrange" 
done 

split(FILENAME, f, /\./)を試してみて、ファイル名から染色体識別子部分を取得するためにf[1]を印刷します。

0

これは私のために働いてしまったものです:

module load SAMtools 

awk '{print $1, $2-5 "-" $2+5}' FILE PATH |\ 
while read chrom range 
do 

    samtools mpileup -f /REFERENCE GENOME\ 
      /${chrom}.COLON BAM FILE\ 
      /${chrom}.MUSCLE BAM FILE\ 
      -r $chrom:$range -o ${chrom}.colon.${range}.pileup 

done

はあなたの助けをありがとう!

関連する問題