☰

タブ区切りバージョン

tag-icon • tag-icon awk sed grep pcregrep
タブ区切りバージョン

以下のように2つの大きなファイルがあります。

ファイル1：

NW_006502347.1 316684
NW_006527876.1 351
NW_006502151.1 27628
NW_006526579.1 232
NW_006525259.1 132
NW_006501641.1 437014
NW_006525259.1 378
NW_006523082.1 215
NW_006522424.1 153
NW_006522101.1 815
NW_006521985.1 505
NW_006521985.1 527
NW_006521722.1 920
NW_006521525.1 73
NW_006521432.1 258
NW_006521302.1 938
NW_006521272.1 585
NW_006521272.1 745
NW_006521038.1 202
NW_006519846.1 1528
NW_006519837.1 10215

ファイル2：

NW_006502347.1  Gnomon  CDS   305319  305340   .  +  0  ID=cds43608  Parent=rna48098  Dbxref=GeneID:102908761,Genbank:XP_006997436.2  Name=XP_006997436.2  gbkey=CDS  gene=LOC102908761  partial=true  product=histone deacetylase 4-like  protein_id=XP_006997436.2
NW_006501037.1  Gnomon  gene  6936    115174   .  -  .  ID=gene0                      Dbxref=GeneID:102922816  Name=Efl1  gbkey=Gene  gene=Efl1  gene_biotype=protein_coding  partial=true  start_range=.,6936
NW_006501037.1  Gnomon  mRNA  6936    115174   .  -  .  ID=rna0      Parent=gene0     Dbxref=GeneID:102922816,Genbank:XM_006970114.2  Name=XM_006970114.2  gbkey=mRNA  gene=Efl1  model_evidence=Supporting evidence includes similarity to: 5 mRNAs%2C 7 Proteins%2C and 99%25 coverage of the annotated genomic feature by RNAseq alignments%2C including 16 samples with support for all annotated introns  partial=true  product=elongation factor like GTPase 1  start_range=.,6936  transcript_id=XM_006970114.2
NW_006501037.1  Gnomon  exon  115095  115174   .  -  .  ID=id1       Parent=rna0      Dbxref=GeneID:102922816,Genbank:XM_006970114.2  gbkey=mRNA  gene=Efl1  partial=true  product=elongation factor like GTPase 1  transcript_id=XM_006970114.2
NW_006501037.1  Gnomon  exon  114246  114355   .  -  .  ID=id2       Parent=rna0      Dbxref=GeneID:102922816,Genbank:XM_006970114.2  gbkey=mRNA  gene=Efl1  partial=true  product=elongation factor like GTPase 1  transcript_id=XM_006970114.2
NW_006502347.1  Gnomon  mRNA  272091  477077   .  +  .  ID=rna48098  Parent=gene26399 Dbxref=GeneID:102908761,Genbank:XM_006997374.2  Name=XM_006997374.2  end_range=477077,.  gbkey=mRNA  gene=LOC102908761  model_evidence=Supporting evidence includes similarity to: 1 mRNA%2C and 90%25 coverage of the annotated genomic feature by RNAseq alignments  partial=true  product=histone deacetylase 4-like  transcript_id=XM_006997374.2
NW_006501037.1  Gnomon  exon  92339   92472    .  -  .  ID=id5       Parent=rna0      Dbxref=GeneID:102922816,Genbank:XM_006970114.2  gbkey=mRNA  gene=Efl1  partial=true  product=elongation factor like GTPase 1  transcript_id=XM_006970114.2
NW_006501037.1  Gnomon  exon  90969   91106    .  -  .  ID=id6       Parent=rna0      Dbxref=GeneID:102922816,Genbank:XM_006970114.2  gbkey=mRNA  gene=Efl1  partial=true  product=elongation factor like GTPase 1  transcript_id=XM_006970114.2
NW_006501037.1  Gnomon  exon  89261   89475    .  -  .  ID=id7       Parent=rna0      Dbxref=GeneID:102922816,Genbank:XM_006970114.2  gbkey=mRNA  gene=Efl1  partial=true  product=elongation factor like GTPase 1  transcript_id=XM_006970114.2
NW_006502151.1  Gnomon  exon  26099   27657    .  -  .  ID=id586652  Parent=rna47002  Dbxref=GeneID:102918658,Genbank:XM_006996663.2  gbkey=mRNA  gene=Rftn1  product=raftlin%2C lipid raft linker 1  transcript_id=XM_006996663.2
NW_006501641.1  Gnomon  mRNA  393496  438556   .  +  .  ID=rna40001  Parent=gene21212 Dbxref=GeneID:102913870,Genbank:XM_015986269.1  Name=XM_015986269.1  gbkey=mRNA  gene=LOC102913870  model_evidence=Supporting evidence includes similarity to: 9 mRNAs%2C 5 Proteins%2C and 81%25 coverage of the annotated genomic feature by RNAseq alignments  product=transmembrane protein 189  transcript_id=XM_015986269.1
NW_006501053.1  Gnomon  exon  5104713 5104872  .  +  .  ID=id45206   Parent=rna3590   Dbxref=GeneID:102916769,Genbank:XR_001580019.1  gbkey=misc_RNA  gene=Mycbpap  product=MYCBP associated protein%2C transcript variant X4  transcript_id=XR_001580019.1
NW_006501053.1  Gnomon  exon  5104959 5105062  .  +  .  ID=id45207   Parent=rna3590   Dbxref=GeneID:102916769,Genbank:XR_001580019.1  gbkey=misc_RNA  gene=Mycbpap  product=MYCBP associated protein%2C transcript variant X4  transcript_id=XR_001580019.1
NW_006501053.1  Gnomon  exon  5105698 5105881  .  +  .  ID=id45208   Parent=rna3590   Dbxref=GeneID:102916769,Genbank:XR_001580019.1  gbkey=misc_RNA  gene=Mycbpap  product=MYCBP associated protein%2C transcript variant X4  transcript_id=XR_001580019.1
NW_006501053.1  Gnomon  exon  5106131 5106246  .  +  .  ID=id45209   Parent=rna3590   Dbxref=GeneID:102916769,Genbank:XR_001580019.1  gbkey=misc_RNA  gene=Mycbpap  product=MYCBP associated protein%2C transcript variant X4  transcript_id=XR_001580019.1

ファイル1の情報を使用して、列13または14から「gene =」の後に続く単語（「Efl1」など）を抽出したいと思います。具体的には、次のことをしたいと思います。

ステップ1)ファイル1の列1のラベル（NW_006527876.1など）をファイル2の列1のラベルと比較し、列1（ファイル2）がファイル1の列1と一致するすべての行を抽出します。

ご覧のとおり、列1（ファイル2）のラベルが重複しているため、ファイル1の各ラベルに複数の一致があります。

ファイル2の列4と5は間隔を表し、列4は間隔の始まり、列5は間隔の終わりです。ファイル 1 の列 2 は、これらの間隔の間の数を表します。

ステップ2) ステップ1で分離された行から、列2（ファイル1）の下の番号がファイル2の列4と5に示す間隔の間にある行を抽出します。

これは私の理解を超えていますが、コマンドの外観は次のとおりです。

awk '{ print $1 }' file 1 |  
awk `$4 *(file2)* < $2 *(file1)*' | awk '$5 *(file2)* > $2 *(file1)*' > output.tsv

出力には、列1の下に固有のラベルを持つ行を含める必要があります。

ステップ3）上記で生成されたoutput.tsvファイルから、列13または14（以下を参照）の下の「gene =」等号の後の単語を抽出して、等号の後の単語のみを含むファイルを作成したいと思います。

最終出力ファイル（この例に基づいています）：

LOC102908761
Rftn1
LOC102913870

答え：

while read -r id pos; do awk -v id="$id" -v pos="$pos" '$1 == id && pos > $4 && pos < $5 { if (gensub(/.*gene=([A-Za-z0-9]*).*/, "\\1", 1) !~ /\s/) print gensub(/.*gene=([A-Za-z0-9]*).*/, "\\1", 1); }' < file2.txt; done <file1.txt > gene_hits.txt

次に、同じ遺伝子の異なる遺伝子座マッピングを維持しながら、重複（複数の間隔の同じ遺伝子座）を削除するには、次のようにします。

perl -ne 'print if ++$k{$_}==1' A_gene_hits.txt > A_genes.txt

ベストアンサー1

区切り文字としてスペースがあるとします。

$ while read -r id pos; do awk -v id="$id" -v pos="$pos" '$1 == id && pos > $4 && pos < $5 { print gensub(/.*gene=([A-Za-z0-9]*).*/, "\\1", 1) }' <file2; done <file1
LOC102908761
Rftn1
LOC102913870

説明する

  • while read -r id pos; do FOO; done <file1：これはfile11行ずつ読み、最初のフィールド（たとえばNW_006502347.1）をシェル変数に入れ、2番目のフィールド$id（たとえば316684）をシェル変数に入れます$pos。その後、FOO各行に対して実行されます。
  • awk -v id="$id" -v pos="$pos" 'BAR' <file2：の各行に対して実行するコマンドをfile1実行します。これにより、一致する部品が検索されます。このスクリプトには、シェルから2つの「外部」変数を知らせる必要があります。つまり、awk変数にはシェル変数と同じ値が割り当てられ、awk変数とシェル変数にも同じ値が割り当てられます。awkBARfile2awkid$idpos$pos
  • $1 == id && pos > $4 && pos < $5：これはスクリプトの「条件付き」部分ですawk。これらの条件が満たされると、次のコマンドが実行されます。ここでは、最初のフィールドが$1現在の行と同じで、のfile24番目と5番目のフィールドの間にあることを確認します。idfile1pos$4$5file2
  • { print gensub(/.*gene=([A-Za-z0-9]*).*/, "\\1", 1) }：上記の条件が満たされると、このコードが実行されます。私たちはそれを最初に変えたいと思いますgensub。gene=英数字文字列の後にランダムな長さが続く文字列を検索します([A-Za-z0-9]*)。英数字の文字列は(括弧で囲まれます。また、文字列全体の前後のすべての文字を)「検索」します。したがって、これは行全体を「検索」し（最初で唯一の）キャプチャグループである次の英数字文字列に置き換えます。最終的に最初の項目を置き換えることを意味しますが、1行に一致が1つしかないと仮定するため、これは意味がありません。.*gene=([A-Za-z0-9]*)"\\1"gene=1gene=

タブ区切りバージョン

一般的に、私はタブで区切られたファイルを使用することを好みます。これにより、特にフィールド9でスペースを区別できます。

while IFS=$'\t' read -r id pos; do awk -F'\t' -v id="$id" -v pos="$pos" '$1 == id && pos > $4 && pos < $5 { print gensub(/.*gene=([A-Za-z0-9]*).*/, "\\1", 1) }' <file2.tsv ; done <file1.tsv

スクリプトの変更は、タブでシェル行とIFS=$'\t'行awkを明示的に分割することです-F'\t'。

おすすめ記事