打印特定数量的核苷酸答案

【问题标题】：Printing a specific number of nucleotides打印特定数量的核苷酸
【发布时间】：2021-08-26 07:21:58
【问题描述】：

我有一个杂合子和纯合子情况的 vcf 统计信息，我想找到与我的 maf 文件匹配的内容。问题是maf文件中的参考字段不同，它排除了替代状态的核苷酸，例如如果您有 ref CAA 并且替代变体是 CAAAAA，则在 maf 文件中您的 ref 将是 AAA。

所以我需要一个代码来更改我的文件中的 ref 字段和 alt 统计信息（可以添加单独的列 ref2 和 alt2）

这是我的文件的 sn-p：

CHR POS ID REF ALT chr11 71579744 rs71049992 A ACAGCAGCTGGACTGGGAGCAGCAGGACCTG（插入案例）

chr11 124880551 rs71859853 CCGGAGT C（删除大小写）

我认为我应该首先计算第 4 列和第 5 列中的核苷酸数。然后如果第 4 列中的数字大于 5（意味着删除），那么在我的 ref2 中，该位置将从下一个与替代核苷酸不同的核苷酸开始。

对于插入，我将更改 alt 位点并跳过参考核苷酸

因此，我想要这个：

CHR 位置 ID REF ALT REF2 ALT2 chr11 71579744 rs71049992 A ACAGCAGCTGGACTGGGAGCAGCAGGACCTG A CAGCAGCTGGACTGGGAGCAGCAGGACCTG

chr11 124880551 rs71859853 CCGGAGT C CGGAGT C

非常感谢您！

【问题讨论】：

你说的是一个maf文件，但是你展示的sn-p看起来不像MAF。
嗨 Armali，是的，这是来自 vcf 统计数据的 sn-p。我想在这个文件和我的 maf 之间找到一个匹配项
对不起，应该是两行 CHR POS ID REF ALT（新行）chr11 71579744 rs71049992 A ACAGCAGCTGGACTGGGAGCAGCAGGACCTG
@Armali。感谢您的友好回复！在更好的示例中，有 8 列，因此问题是以某种方式仅从 4 列和 5 列中获取不匹配的核苷酸，并将它们与 maf 文件中的 ref 列匹配
@Armali 是的，完全正确！这将是删除的情况。对于插入，我的 ALT2 列应该包含不匹配的核苷酸。我想我应该首先计算我的 vcf 统计文件中哪一列（4 或 5）有更多的核苷酸？

标签： reference vcf-vcard maf

【解决方案1】：

我想我应该先数一下第 4 列和第 5 列的核苷酸数……

使用awk，您可以使用length 函数计算核苷酸的数量：

awk 'NR==1 {print $0" REF2 ALT2"}       # assuming first line has column headers
     NR>1  {if (length($4)<length($5)) print $0, $4, gensub($4, "", 1, $5)
            else                       print $0, gensub($5, "", 1, $4), $5}' file

【讨论】：

非常感谢！目前试图将代码理解为无法正常工作。 awk 'NR==1 {print $0" REF2 ALT2"} NR>1 {if (length($4)1 {if (length($4)
如果else不在行首，则;必须放在前面。
完美运行！非常感谢！！
只是好奇，请你解释一下这部分 gensub($4, "", 1, $5)
gensub($4, "", 1, $5) 用"" 替换$4（第4 列）中$4（第4 列）的匹配项""（什么都没有，即删除匹配项）；第 5 列未更改，但函数返回修改后的字符串。