http://genome.sph.umich.edu/w/index.php?action=history&feed=atom&title=830_-_MarkovModel%3A%3AImpute%28%29 2026-09-27T23:05:08Z Revision history for this page on the wiki MediaWiki 1.43.1 http://genome.sph.umich.edu/w/index.php?title=830_-_MarkovModel::Impute()&diff=8734&oldid=prev 2013-09-25T18:50:19Z

<p>Created page with &quot;&lt;source lang=&quot;cpp&quot;&gt; void MarkovModel::Impute(char * major, char * observed, float * probs, char ** haplotypes, float ** freqs, int position) { d...&quot;</p> <p><b>New page</b></p><div>&lt;source lang=&quot;cpp&quot;&gt;<br /> void MarkovModel::Impute(char * major, char * observed, float * probs,<br /> char ** haplotypes, float ** freqs, int position)<br /> {<br /> double P[5] = {0.0, 0.0, 0.0, 0.0, 0.0};<br /> <br /> for (int i = 0; i &lt; states; i++)<br /> P[haplotypes[i][position]] += probs[i];<br /> <br /> double ptotal = P[1] + P[2] + P[3] + P[4];<br /> double pmajor = P[major[position]];<br /> <br /> int mle = 1;<br /> for (int i = 2; i &lt;= 4; i++)<br /> if (P[i] &gt;= P[mle])<br /> mle = i;<br /> <br /> char labels[] = {0, &#039;a&#039;, &#039;c&#039;, &#039;g&#039;, &#039;t&#039;};<br /> <br /> imputedDose[position] += imputedHap[position] = (pmajor / ptotal);<br /> imputedAlleles[position] = labels[mle];<br /> <br /> // What would we have imputed if we didn&#039;t observe a genotype at the current position?<br /> // This assumes that freqs[i][position] stores the frequency of allele i at position (for i &gt; 0)<br /> // and is equal to 1.0 when the current genotype is not observed (i = 0)<br /> <br /> double fmatch = 1.0 / (1. - E[position] + E[position] * freqs[observed[position]][position]);<br /> double fmismatch = 1.0 / (E[position] * freqs[observed[position]][position]);<br /> <br /> for (int i = 1; i &lt;= 4; i++)<br /> if (observed[position] == i)<br /> P[i] *= fmatch;<br /> else<br /> P[i] *= fmismatch;<br /> <br /> ptotal = P[1] + P[2] + P[3] + P[4];<br /> pmajor = P[major[position]];<br /> <br /> leaveOneOut[position] = pmajor / ptotal;<br /> }<br /> &lt;/source&gt;</div>

Goncalo