確立された:File_A と重なる File_B から名前を抽出します。
file_bの2番目の列が4列とfile_b5の間にある場合は、file_bの最初の列をfile_a(通常「Name =」の後の列10)と一致させて遺伝子名を抽出したいと思います。 1列が一致する必要があるため、行(file_b)ごとに1つの遺伝子のみを取得します。しかし、理論的には、同じ遺伝子と一致する複数の隣接する行(column_b)を持つことができます(例:file_bの2行目MT 4065
)。
現在のコードにはいくつかの問題があります。
(1)以下の設定方法でfile_bの最後の行が出力から欠落していますが、そのgroupVII 17978350
行()がリストに表示されている場合は状況が異なります。設定どおりに機能したいです。
(2)名前に特殊文字(コロン、ハイフンなど)が含まれている場合、名前は切り捨てられます。等号の後にフルネームを追加したい。
(3)最初の2列が項目で、3列目が遺伝子ヒットになるように、file_bの項目/行を出力の遺伝子ヒットと一致させたいと思います。
file_a.tsv
MT insdc gene 2851 3825 . + . ID=gene:ENSGACG00000020925 Name=mt-nd1 biotype=protein_coding description=NADH dehydrogenase 1%2C mitochondrial [Source:ZFIN%3BAcc:ZDB-GENE-011205-7] gene_id=ENSGACG00000020925 logic_name=mt_genbank_import version=1
MT insdc gene 4036 5082 . + . ID=gene:ENSGACG00000020929 Name=mt-nd2 biotype=protein_coding description=NADH dehydrogenase 2%2C mitochondrial [Source:ZFIN%3BAcc:ZDB-GENE-011205-8] gene_id=ENSGACG00000020929 logic_name=mt_genbank_import version=1
groupIII ensembl gene 7332324 7334769 . - . ID=gene:ENSGACG00000015265 Name=si:dkeyp-68b7.10 biotype=protein_coding description=si:dkeyp-68b7.10 [Source:ZFIN%3BAcc:ZDB-GENE-070912-667] gene_id=ENSGACG00000015265 logic_name=ensembl version=1
groupIV ensembl gene 1368026 1374881 . + . ID=gene:ENSGACG00000016447 Name=hnrnpa0b biotype=protein_coding description=heterogeneous nuclear ribonucleoprotein A0b [Source:ZFIN%3BAcc:ZDB-GENE-030131-6154] gene_id=ENSGACG00000016447 logic_name=ensembl version=1
groupIV ensembl gene 5347339 5349041 . - . ID=gene:ENSGACG00000017010 Name=zgc:153018 biotype=protein_coding description=zgc:153018 [Source:ZFIN%3BAcc:ZDB-GENE-060929-752] gene_id=ENSGACG00000017010 logic_name=ensembl version=1
groupV ensembl gene 120615 125489 . + . ID=gene:ENSGACG00000002103 Name=zdhhc6 biotype=protein_coding description=zinc finger%2C DHHC-type containing 6 [Source:ZFIN%3BAcc:ZDB-GENE-030131-3189] gene_id=ENSGACG00000002103 logic_name=ensembl version=1
groupVI ensembl gene 11230354 11232784 . + . ID=gene:ENSGACG00000009527 Name=bnip4 biotype=protein_coding description=BCL2 interacting protein 4 [Source:ZFIN%3BAcc:ZDB-GENE-051113-212] gene_id=ENSGACG00000009527 logic_name=ensembl version=1
groupVII ensembl gene 2271611 2277214 . + . ID=gene:ENSGACG00000019012 Name=sf3b2 biotype=protein_coding description=splicing factor 3b%2C subunit 2 [Source:ZFIN%3BAcc:ZDB-GENE-070928-1] gene_id=ENSGACG00000019012 logic_name=ensembl version=2
groupVII ensembl gene 15815857 15824549 . + . ID=gene:ENSGACG00000020296 Name=mpp1 biotype=protein_coding description=membrane protein%2C palmitoylated 1 [Source:ZFIN%3BAcc:ZDB-GENE-031113-4] gene_id=ENSGACG00000020296 logic_name=ensembl version=1
groupVII ensembl gene 17978322 17982388 . + . ID=gene:ENSGACG00000020399 Name=si:ch211-284e13.4 biotype=protein_coding description=si:ch211-284e13.4 [Source:ZFIN%3BAcc:ZDB-GENE-060526-161] gene_id=ENSGACG00000020399 logic_name=ensembl version=1
ファイル_b.tsv
MT 4050
groupIII 7332350
groupIV 5347350
groupVI 11230375
groupVII 17978350
パスワード:
while read -r id pos; do awk -v id="$id" -v pos="$pos" '$1 == id && pos > $4 && pos < $5 { if (gensub(/.*Name=([A-Za-z0-9]*).*/, "\\1", 1) !~ /\s/) print gensub(/.*Name=([A-Za-z0-9]*).*/, "\\1", 1); }' <file_a.tsv; done < file_b.tsv > output.tsv
出力.tsv
mt
si
zgc
bnip4
期待される出力
MT 4050 mt-nd2
groupIII 7332350 si:dkeyp-68b7.10
groupIV 5347350 zgc:153018
groupVI 11230375 bnip4
groupVII 17978350 si:ch211-284e13.4
答え1
# save this as script.awk or whatevernameyouwant.awk
function within_range(val, lower, upper, proximity) {
# you can specify the "proximity" as required
return val > lower - proximity && val < upper + proximity
}
BEGIN {
OFS="\t"
}
$1 == id && within_range(pos, $4, $5, 100) {
name = gensub(/.*Name=([^\t]*).*/, "\\1", 1)
if (name ~ /[^[:space:]]+/)
print id, pos, name
}
その後実行
while read -r id pos
do
awk -v id=$id -v pos=$pos -f script.awk file_a.tsv
done < file_b.tsv > output.tsv
お願いします.tsv
ファイルを処理する前に、ファイルのフィールドをタブで区切る必要があります。私の結果:
MT 4050 mt-nd2
groupIII 7332350 si:dkeyp-68b7.10
groupIV 5347350 zgc:153018
groupVI 11230375 bnip4
groupVII 17978350 si:ch211-284e13.4
IDの場合、遺伝子ヒットは不可能でなければMT
なりません。mt-nd2
mt-nd1
データ処理にはまだPythonを使用することをお勧めします。
答え2
予想される表示出力が私の意見に一致しないようです(2行 - > 1行と3行)。これがオタだったら、次のことを試すことができますか?
awk 'FNR==NR{a[$1]=$2;next} ($1 in a) && (a[$1]>=$4 && a[$1]<=$5){sub("Name=","",$10);print $1,a[$1],$10}' b.tsv a.tsv > output.tsv