根据后验平均基因型计算等位基因频率答案

【问题标题】：Calculate allele frequencies from posterior mean genotypes根据后验平均基因型计算等位基因频率
【发布时间】：2021-05-08 04:47:48
【问题描述】：

在以下文件（input_file.txt）中：

col1, col2, col3, NA, NA, NA, NA, 0, 0.1, 0, 0.002, 1, 0, NA, NA, NA, 1, 1.9, 0, 0.1, 2, 1, 0, 0.8

后平均基因型 (PMG) 为 0-2，或 NA（前 3 列不包含 PMG）。

因此，input_file.txt 包含 21 个人的 PMG 信息（=字段数 (NF) -3）。

要根据 input_file.txt 中包含的信息计算等位基因频率，我们首先需要进行一些预处理并将一些值保存在内存中。

如上所述，PMG 的范围为 0-2，但这些值中的每一个都代表仅有的三种基因型之一：

0 到
0.5 到
1.5 到 2 = 纯合替代

所以我们接下来需要将 PMG 转换为 input_file.txt 中的基因型

类似：

for NF >3 && column value !=NA 
if "column value" = 0 to <0.5, replace value with 0 (= homozygous reference)
if "column value" = 0.5 to <1.5, replace value with 1 (= heterozygous)
if "column value" = 1.5 to 2, replace value with 2 (= homozygous alternate)

中间输出：

col1, col2, col3, NA, NA, NA, NA, 0, 0, 0, 0, 1, 0, NA, NA, NA, 1, 2, 0, 0, 2, 1, 0, 1

所以我们为 21 个人中的 14 个人提供了非北美 PMG。

我可以使用以下方法确定计数：

awk -F, ' NF>3 {print (split($0,a,"NA")-1) }' input_file.txt > NA_count
awk -F, ' NF>3 {print (split($0,a,"0")-1) }' input_file.txt > ref_count
awk -F, ' NF>3 {print (split($0,a,"1")-1) }' input_file.txt > het_count
awk -F, ' NF>3 {print (split($0,a,"2")-1) }' input_file.txt > alt_count

要从此数据中获取等位基因频率： (ref_count + het_count + 2alt_count)/2(ref_count + het_count + alt_count)

等同于： 0+0+0+0+1+0+1+2+0+0+2+1+0+1/2(14) = 8/28 = 0.286 = 等位基因频率

期望的输出：

"col1", 8/28, 0.286

请注意，输入文件可以有可变数量的列，其中会有可变数量的 0-2 和 NA 值。

上面显示的中间输出本身不是必需的（只需要所需的输出）

awk 代码首选但不是必需的

【问题讨论】：

标签： arrays awk bioinformatics

【解决方案1】：

awk 来救援！

$ awk -F', *' -v OFS=', ' '
       {for(i=4;i<=NF;i++) if($i!="NA") a[$i=int($i+0.5)]++}
   1;
   END {print "\""$1"\"", (n=a[1]+2*a[2])"/"(d=2*(a[0]+a[1]+a[2])),n/d}' file

col1, col2, col3, NA, NA, NA, NA, 0, 0, 0, 0, 1, 0, NA, NA, NA, 1, 2, 0, 0, 2, 1, 0, 1
"col1", 8/28, 0.285714

要过滤中间结果，请从代码中删除1;。

【讨论】：