18#include "GenotypeLists.h"
24GenotypeList::GenotypeList()
29bool GenotypeList::EliminateGenotypes(
Pedigree & ped,
Family * family,
int marker)
35 InitializeList(list, ped, family, marker);
40#ifdef DEBUG_ELIMINATOR
41 Print(list, ped, family, marker);
44 while (PairwiseCheck(list, ped, family) || FamilyCheck(list, ped, family))
45#ifdef DEBUG_ELIMINATOR
46 Print(list, ped, family, marker)
50 for (
int i = 0; i < family->count; i++)
51 if (!list[i].ignore && list[i].allele1.Length() == 0)
53 printf(
"%s - Family %s has a subtle genotype inconsistency\n",
54 (
const char *) ped.markerNames[marker], (
const char *) family->famid);
66 for (
int i = family->count - 1; i >= 0; i--)
68 Person & person = ped[family->path[i]];
69 int id = person.traverse;
70 bool maleX = person.sex == SEX_MALE && ped.chromosomeX;
72#ifdef DEBUG_ELIMINATOR
73 printf(
"Initializing genotype list for %s ...\n", (
const char *) person.pid);
77 if (person.markers[marker].isKnown())
80 list[id].Dimension(1);
81 list[id].SetGenotype(0, person.markers[marker][0], person.markers[marker][1]);
82 list[id].alleles.Clear();
83 list[id].alleles.Push(person.markers[marker][0]);
84 list[id].alleles.PushIfNew(person.markers[marker][1]);
85 list[id].ignore =
false;
88 if (maleX && person.markers[marker].isHeterozygous())
89 list[id].Dimension(0);
91 else if (list[
id].alleles.Length())
92 if (person.sex == SEX_MALE && ped.chromosomeX)
95 list[id].Dimension(list[
id].alleles.Length() + 1);
97 for (
int i = 0, out = 0; i < list[id].alleles.Length(); i++)
98 list[
id].SetGenotype(out++, list[
id].alleles[i], list[
id].alleles[i]);
99 list[id].SetGenotype(list[
id].alleles.Length(), -1, -1);
101 list[id].ignore =
false;
106 int count = list[id].alleles.Length() * (list[id].alleles.Length() + 3) / 2 + 1;
108 list[id].Dimension(count);
110 for (
int i = 0, out = 0; i < list[id].alleles.Length(); i++)
113 for (
int j = 0; j <= i; j++)
114 list[
id].SetGenotype(out++, list[
id].alleles[i], list[
id].alleles[j]);
117 list[id].SetGenotype(out++, list[
id].alleles[i], -1);
121 list[id].SetGenotype(count - 1, -1, -1);
123 list[id].ignore =
false;
126 list[id].ignore =
true;
129 if (i < family->founders)
continue;
132 int fatid = person.father->traverse;
133 int motid = person.mother->traverse;
135 for (
int i = 0; i < list[id].alleles.Length(); i++)
137 list[motid].alleles.PushIfNew(list[
id].alleles[i]);
138 if (!maleX) list[fatid].alleles.PushIfNew(list[
id].alleles[i]);
145#ifdef DEBUG_ELIMINATOR
146 printf(
"Checking Relative Pairs ...\n");
149 bool changed =
false;
151 for (
int i = family->count - 1; i >= family->founders; i--)
153 Person & person = ped[family->path[i]];
155 int id = person.traverse;
156 int fatid = person.father->traverse;
157 int motid = person.mother->traverse;
159 bool maleX = person.sex == SEX_MALE && ped.chromosomeX;
161 if (list[
id].ignore)
continue;
164 for (
int i = 0; i < list[id].allele1.Length(); i++)
166 int al1 = list[id].allele1[i];
167 int al2 = list[id].allele2[i];
170 if ((maleX && !list[motid].Matches(al1) && al1 != -1) ||
171 (!maleX && !(al1 == -1 && al2 == -1) &&
172 !(list[fatid].Matches(al1) && (al2 == -1 || list[motid].Matches(al2))) &&
173 !((al2 == -1 || list[fatid].Matches(al2)) && list[motid].Matches(al1))))
175 list[id].Delete(i--);
182 if (list[
id].Matches(-1))
186 for (
int i = 0; i < list[motid].allele1.Length(); i++)
188 int al1 = list[motid].allele1[i];
189 int al2 = list[motid].allele2[i];
192 if (!list[
id].Matches(al1) &&
193 !list[id].Matches(al2))
195 list[motid].Delete(i--);
204 for (
int i = 0; i < list[fatid].allele1.Length(); i++)
206 int al1 = list[fatid].allele1[i];
207 int al2 = list[fatid].allele2[i];
210 if (!list[
id].Matches(al1) &&
211 !list[id].Matches(al2))
213 list[fatid].Delete(i--);
218#ifdef DEBUG_ELIMINATOR
219 printf(
"Done checking individual %s\n", (
const char *) person.pid);
220 Print(list, ped, family, 0);
230#ifdef DEBUG_ELIMINATOR
231 printf(
"Checking Nuclear Families ...\n");
234 bool changed =
false;
236 for (
int i = family->count - 1; i >= family->founders; i--)
238 Person & person = ped[family->path[i]];
240 int fatid = person.father->traverse;
241 int motid = person.mother->traverse;
244 if (person.sibs[0] != &person || list[fatid].ignore || list[motid].ignore)
247#ifdef DEBUG_ELIMINATOR
248 printf(
"Checking Sibship with %s ...\n", (
const char *) person.pid);
252 list[fatid].checked = 0;
253 list[motid].checked = 0;
255 for (
int i = 0; i < person.sibCount; i++)
256 list[person.sibs[i]->traverse].checked = 0;
259 changed |= TrimParent(list, person, fatid, motid);
262 changed |= TrimParent(list, person, motid, fatid);
265 for (
int i = 0; i < person.sibCount; i++)
267 int sibid = person.sibs[i]->traverse;
268 bool maleX = person.sibs[i]->sex == SEX_MALE && ped.chromosomeX;
273 for (
int j = list[sibid].checked; j < list[sibid].allele1.Length(); j++)
274 changed |= Cleanup(list, person, motid, fatid, sibid, j);
277#ifdef DEBUG_ELIMINATOR
285bool GenotypeList::Matches(
int genotype,
int allele)
287 return allele1[genotype] == allele || allele2[genotype] == allele;
290bool GenotypeList::Matches(
int allele)
292 return allele1.Find(allele) != -1 || allele2.Find(allele) != -1;
295int GenotypeList::SaveGenotype(
int genotype)
297 if (checked > genotype)
300 if (checked != genotype)
302 allele1.Swap(genotype, checked);
303 allele2.Swap(genotype, checked);
309bool GenotypeList::CheckTrio(
GenotypeList * list,
int fatid,
int motid,
int child,
313 return (list[fatid].Matches(i, list[child].allele1[k]) &&
314 (list[motid].Matches(j, list[child].allele2[k]) || list[child].allele2[k] == -1)) ||
315 ((list[fatid].Matches(i, list[child].allele2[k]) || list[child].allele2[k] == -1) &&
316 list[motid].Matches(j, list[child].allele1[k])) ||
317 (list[child].allele1[k] == -1 && list[child].allele2[k] == -1);
320void GenotypeList::Dimension(
int genotypes)
322 allele1.Dimension(genotypes);
323 allele2.Dimension(genotypes);
326void GenotypeList::SetGenotype(
int genotype,
int al1,
int al2)
328 allele1[genotype] = al1;
329 allele2[genotype] = al2;
332void GenotypeList::Delete(
int genotype)
334 allele1.Delete(genotype);
335 allele2.Delete(genotype);
338bool GenotypeList::TrimParent(
GenotypeList * list,
Person & person,
int motid,
int fatid)
340 bool trimmed =
false;
342 while (list[motid].checked < list[motid].allele1.Length())
344 int current = list[motid].allele1.Length() - 1;
348 for (
int i = list[fatid].allele1.Length() - 1; i >= 0; i--)
353 for (
int j = 0; j < person.sibCount; j++)
355 int sibid = person.sibs[j]->traverse;
356 int maleX = person.sibs[j]->sex == SEX_MALE && person.chromosomeX;
360 if (list[sibid].ignore || maleX)
366 for (
int k = list[sibid].allele1.Length() - 1; k >= 0; k--)
367 if (CheckTrio(list, motid, fatid, sibid, current, i, k))
373 if (matches != j + 1)
378 if (matches == person.sibCount)
380 for (
int j = 0; j < person.sibCount; j++)
382 int sibid = person.sibs[j]->traverse;
384 for (
int k = list[sibid].checked; k < list[sibid].allele1.Length(); k++)
385 if (CheckTrio(list, motid, fatid, sibid, current, i, k))
386 list[sibid].SaveGenotype(k);
389 list[motid].SaveGenotype(current);
390 list[fatid].SaveGenotype(i);
400 list[motid].Delete(current);
408bool GenotypeList::Cleanup(
GenotypeList * list,
Person & person,
int motid,
int fatid,
int child,
int geno)
410 for (
int current = 0; current < list[motid].allele1.Length(); current++)
411 for (
int i = list[fatid].allele1.Length() - 1; i >= 0; i--)
412 if (CheckTrio(list, motid, fatid, child, current, i, geno))
417 for (
int j = 0; j < person.sibCount; j++)
419 int sibid = person.sibs[j]->traverse;
420 int maleX = person.sibs[j]->sex == SEX_MALE && person.chromosomeX;
424 if (list[sibid].ignore || maleX)
430 for (
int k = list[sibid].allele1.Length() - 1; k >= 0; k--)
431 if (CheckTrio(list, motid, fatid, sibid, current, i, k))
437 if (matches != j + 1)
442 if (matches == person.sibCount)
443 for (
int j = 0; j < person.sibCount; j++)
445 int sibid = person.sibs[j]->traverse;
447 for (
int k = list[sibid].checked; k < list[sibid].allele1.Length(); k++)
448 if (CheckTrio(list, motid, fatid, sibid, current, i, k))
449 list[sibid].SaveGenotype(k);
455 list[child].Delete(geno);
462 MarkerInfo * info = ped.GetMarkerInfo(marker);
464 for (
int i = 0; i < family->count; i++)
466 printf(
"%s - ", (
const char *) ped[family->path[i]].pid);
468 for (
int j = 0; j < list[i].allele1.Length(); j++)
470 if (list[i].allele1[j] == -1)
473 printf(
"%s/", (
const char *) info->GetAlleleLabel(list[i].allele1[j]));
475 if (list[i].allele2[j] == -1)
478 printf(
"%s ", (
const char *) info->GetAlleleLabel(list[i].allele2[j]));