<p>Beside the correct question by @Kamil, this sounds like a job for <a href="https://www.gnu.org/software/datamash/" rel="nofollow noreferrer">datamash</a>.</p>

<p>First let's get the number of columns in the file</p>

<pre><code>$ head -n1 genotypes.txt| awk -v FS="\t" '{print NF}'
</code></pre>

<p>Now lets calculate the average per column with <code>datamash</code></p>

<pre><code>$ datamash -H mean 2-&lt;n&gt; &lt; genotypes.txt &gt; mean.txt
</code></pre>

<p>Replace <code>&lt;n&gt;</code> with the number of columns we received with command above.</p>

<p>For your example the output will look like this:</p>

<pre><code>mean(SNP1)  mean(SNP2)  mean(SNP3)  mean(SNP4)
0,5 0   1,5 1,5
</code></pre>

<p>One can remove the <code>mean</code> with <code>sed</code></p>

<pre><code>$ sed -i -e 's/mean(//g' -e 's/)//g' mean.txt
</code></pre>

CC BY-SA 4.0
2019-02-08T14:21:01.823Z
Sorry, this is my fault for not providing good enough example data - there are many pairs of rows in the file - will edit the example data to make this clear.
Joni
CC BY-SA 4.0
2019-02-08T12:05:13.173Z
2019-02-08T12:05:13.173Z
finswimmer
<p>These are unusual parents because they all appear haploid. Anyway, if I was doing this I would calculate the allele frequency per progeny such as makefreq in adegenet in R-cran (just my personal thoughts). </p>

<p><strong>Problem</strong> The issue is that you are using a continuous distribution 1,2,3 etc... for catagorical data and will result in confusing homozygotes with heterozygotes. </p>

<p><strong><em>For example</em></strong></p>

<ul>
<li>If you have allele 3 (father) and allele 1 (mother) the "progeny" wil
be "allele" 2</li>
<li>If you have allele 2 (father) and allele 2 (mother) the "progeny"
will be "allele" 2</li>
</ul>

<p>However, one is homozygous and the other is heterozygous and the allelic background for each is quite differnt. </p>

<p>The alternative is to use letters rather than numbers. You can always use "OneHotKeyEncoder" to convert back from AA, aa etc... to numerical catagories for analysis later (its a standard feature in Python).</p>

<p>If you used letters rather and then converted to numbers every "progeny" would have its own numerical catagory (but you would lose homozygous vs. heterozygous classification). It is standard practice to switch from letters to numbers and the "encoder" will ensure no mistakes are made. The other way around this is to designate the variable as a string, rather than as an integer... this might be more to your liking, but you are stuck in base 10.</p>

<p>With regards the throughput, 100000s 1000 individuals, a 100 million SNPs is doable via R, the operations are simple.</p>

<p>Just to clarify: it would be easy to script to a solution (even from base code), but I'm not sure about how to interpret the output. All I'm suggesting is you have another think about the goal of this analysis. </p>

CC BY-SA 4.0
2019-02-11T17:01:50.270Z
Again, my lack of clarity - 0,1,2 relates to B allele counts of bivariate SNPs (A/B) - these are not categorical values
Joni
CC BY-SA 4.0
2019-02-11T21:20:38.603Z
In that case you appear to be violating Mendelian inheritance. Children can't receive 0.5 allele, which is what you are implying
M__
CC BY-SA 4.0
2019-02-12T22:33:28.657Z
No - they either receive 1, or 0. But on average, they receive 0.5.
Joni
CC BY-SA 4.0
2019-02-08T15:23:41.230Z
2019-02-08T15:31:53.057Z
2019-02-08T15:31:53.057Z
M__
<p>I am interested in creating an averaged value for genotype across two parents:</p>

<p>E.g. Input (three unrelated couples, 1, 2 and 3):</p>

<pre><code>Individual    SNP1    SNP2    SNP3    SNP4 ...
Father_1      1       0       1       2    ...
Mother_1      0       0       2       1    ...
Father_2      1       0       1       1    ...
Mother_2      2       1       0       1    ...
Father_3      1       2       0       1    ...
Mother_3      1       .       0       2    ...
</code></pre>

<p>E.g. Output (the midparental genotypes for couple 1, 2 and 3):</p>

<pre><code>Individual    SNP1    SNP2    SNP3    SNP4 ...
Midparent_1   0.5     0       1.5     1.5  ...
Midparent_2   1.5     0.5     0.5     1    ...
Midparent_3   1       .       0       1.5  ...
</code></pre>

<p>I can think of a way to do this using e.g. a PLINK .raw file and manipulating that in R, but am concerned that would become memory- and storage-heavy when assessing 100000s of SNPs and 1000s of pairs.</p>

<p>Does anyone know of an existing tool that can do this?</p>

<p>Thanks</p>

CC BY-SA 4.0
2019-02-08T10:33:33.520Z
Hi there, maybe it's just my ignorance, but what "average parental genotype" should represent? Two homozygotes (0 2) and two heterozygotes (1 1) have the same average (1) but produce quite different genotypes (0 1 2 in one case 1 in the other). Performance wise, would guess that a small python script could do this quite fast. 100,000 SNPs is not that much.
Kamil S Jaron
CC BY-SA 4.0
2019-02-08T10:34:25.220Z
Welcome to the site Joni. Could you show us how would you do it using PLINK and R? This way we could improve the code or suggest faster/more efficient alternatives knowing the base line. To do so, remember you can [edit] your question.
llrs
CC BY-SA 4.0
2019-02-08T14:13:35.977Z
Question is also posted on [biostars](https://www.biostars.org/p/362754/).
finswimmer
CC BY-SA 4.0
2019-02-08T14:26:41.967Z
Example pseudocode of how this could be achieved in PLINK and R. 1) Split PLINK binary file to couples and offspring. 2) Convert PLINK couple binary file to a .raw file with --recode A. 3) Sort by FID (assuming FID identifies couples). 4) Create new PLINK midparent.raw file in R by averaging each SNP across the two lines. 5) Create additional files needed to convert between PLINK .raw and PLINK binary. 6) Create new PLINK binary with midparent.raw. 7) Merge with offspring binary.
Joni
<snp><plink><data-management>
Calculate average parental genotype
# Calculate average parental genotype

(Asked by: username_0 on 2019-02-08)

I am interested in creating an averaged value for genotype across two parents:
E.g. Input (three unrelated couples, 1, 2 and 3):
Individual    SNP1    SNP2    SNP3    SNP4 ...
Father_1      1       0       1       2    ...
Mother_1      0       0       2       1    ...
Father_2      1       0       1       1    ...
Mother_2      2       1       0       1    ...
Father_3      1       2       0       1    ...
Mother_3      1       .       0       2    ...

E.g. Output (the midparental genotypes for couple 1, 2 and 3):
Individual    SNP1    SNP2    SNP3    SNP4 ...
Midparent_1   0.5     0       1.5     1.5  ...
Midparent_2   1.5     0.5     0.5     1    ...
Midparent_3   1       .       0       1.5  ...

I can think of a way to do this using e.g. a PLINK .raw file and manipulating that in R, but am concerned that would become memory- and storage-heavy when assessing 100000s of SNPs and 1000s of pairs.
Does anyone know of an existing tool that can do this?
Thanks

**Tags:** snp, plink, data-management
> **Comments:**
> **username_1**: Hi there, maybe it's just my ignorance, but what "average parental genotype" should represent? Two homozygotes (0 2) and two heterozygotes (1 1) have the same average (1) but produce quite different genotypes (0 1 2 in one case 1 in the other). Performance wise, would guess that a small python script could do this quite fast. 100,000 SNPs is not that much.  
> **username_2**: Welcome to the site Joni. Could you show us how would you do it using PLINK and R? This way we could improve the code or suggest faster/more efficient alternatives knowing the base line. To do so, remember you can [edit] your question.  
> **username_3**: Question is also posted on [biostars](https://www.biostars.org/p/362754/).  
> **username_0**: Example pseudocode of how this could be achieved in PLINK and R. 1) Split PLINK binary file to couples and offspring. 2) Convert PLINK couple binary file to a .raw file with --recode A. 3) Sort by FID (assuming FID identifies couples). 4) Create new PLINK midparent.raw file in R by averaging each SNP across the two lines. 5) Create additional files needed to convert between PLINK .raw and PLINK binary. 6) Create new PLINK binary with midparent.raw. 7) Merge with offspring binary.  

---

## 2 answers

### Answer by username_3 on 2019-02-08. Score: 2

Beside the correct question by @Kamil, this sounds like a job for datamash.
First let's get the number of columns in the file
$ head -n1 genotypes.txt| awk -v FS="\t" '{print NF}'

Now lets calculate the average per column with datamash
$ datamash -H mean 2- < genotypes.txt > mean.txt

Replace  with the number of columns we received with command above.
For your example the output will look like this:
mean(SNP1)  mean(SNP2)  mean(SNP3)  mean(SNP4)
0,5 0   1,5 1,5

One can remove the mean with sed
$ sed -i -e 's/mean(//g' -e 's/)//g' mean.txt

> **Comments:**
> **username_0**: Sorry, this is my fault for not providing good enough example data - there are many pairs of rows in the file - will edit the example data to make this clear.  

---

### Answer by username_4 on 2019-02-08. Score: 1

These are unusual parents because they all appear haploid. Anyway, if I was doing this I would calculate the allele frequency per progeny such as makefreq in adegenet in R-cran (just my personal thoughts).

**Problem** The issue is that you are using a continuous distribution 1,2,3 etc... for catagorical data and will result in confusing homozygotes with heterozygotes.

***For example***

* If you have allele 3 (father) and allele 1 (mother) the "progeny" wil
  be "allele" 2
* If you have allele 2 (father) and allele 2 (mother) the "progeny"
  will be "allele" 2

However, one is homozygous and the other is heterozygous and the allelic background for each is quite differnt.

The alternative is to use letters rather than numbers. You can always use "OneHotKeyEncoder" to convert back from AA, aa etc... to numerical catagories for analysis later (its a standard feature in Python).

If you used letters rather and then converted to numbers every "progeny" would have its own numerical catagory (but you would lose homozygous vs. heterozygous classification). It is standard practice to switch from letters to numbers and the "encoder" will ensure no mistakes are made. The other way around this is to designate the variable as a string, rather than as an integer... this might be more to your liking, but you are stuck in base 10.

With regards the throughput, 100000s 1000 individuals, a 100 million SNPs is doable via R, the operations are simple.

Just to clarify: it would be easy to script to a solution (even from base code), but I'm not sure about how to interpret the output. All I'm suggesting is you have another think about the goal of this analysis.

> **Comments:**
> **username_0**: Again, my lack of clarity - 0,1,2 relates to B allele counts of bivariate SNPs (A/B) - these are not categorical values  
> **username_4**: In that case you appear to be violating Mendelian inheritance. Children can't receive 0.5 allele, which is what you are implying  
> **username_0**: No - they either receive 1, or 0. But on average, they receive 0.5.  

---

