Loading...
Searching...
No Matches
PedigreeTrim.cpp
/*
* Copyright (C) 2010 Regents of the University of Michigan
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#include "Pedigree.h"
void Pedigree::ShowTrimHeader(bool & flag)
{
if (flag)
{
printf("Trimming uninformative individuals...\n");
flag = false;
}
}
void Pedigree::Trim(bool quiet, int * informative)
{
int newCount = 0;
Person ** newPersons = new Person * [count];
// This function applies the following filters to reduce complexity
// of pedigree
//
// RULE 1: Remove all pedigrees no genotype or phenotype data
// RULE 2: Remove leaf individuals with no data
// RULE 3: Remove founder couples with <2 offspring and no data
bool showHeader = true;
IntArray discardable, offspring, mates, haveData;
for (int f = 0; f < familyCount; f++)
{
Family * fam = families[f];
// Cache for storing indicators about whether each family member is
// informative
haveData.Dimension(fam->count);
// Check that some data is available in the family
int hasData = false;
for (int i = fam->first; i <= fam->last; i++)
if (informative == NULL)
hasData |= haveData[persons[i]->traverse] = persons[i]->haveData();
else
hasData |= haveData[persons[i]->traverse] = informative[i];
if (!hasData)
{
if (!quiet)
{
ShowTrimHeader(showHeader);
printf(" Removing family %s: No data\n", (const char *) fam->famid);
}
for (int i = fam->first; i <= fam->last; i++)
delete persons[i];
continue;
}
// Assume that we need everyone in the family
discardable.Dimension(fam->count);
discardable.Set(0);
bool trimming = true;
while (trimming)
{
trimming = false;
// Tally the number of offspring for each individual
offspring.Dimension(fam->count);
offspring.Zero();
// Tally the number of mates for each individual
mates.Dimension(fam->count);
mates.Set(-1);
// In the first round, we count the number of offspring
// for each individual in the current trimmed version of the
// pedigree
for (int i = fam->count - 1; i >= fam->founders; i--)
{
if (discardable[i]) continue;
Person & p = *(persons[fam->path[i]]);
if (discardable[p.father->traverse])
continue;
if (offspring[i] == 0 && !haveData[p.traverse])
{
trimming = true;
discardable[i] = true;
continue;
}
int father = p.father->traverse;
int mother = p.mother->traverse;
if (mates[father] == -1 && mates[mother] == -1)
{
mates[father] = mother,
mates[mother] = father;
}
else if (mates[father] != mother)
{
if (mates[father] >= 0)
mates[mates[father]] = -2;
if (mates[mother] >= 0)
mates[mates[mother]] = -2;
mates[mother] = -2;
mates[father] = -2;
}
offspring[father]++;
offspring[mother]++;
}
// In the second pass, we remove individuals with no
// data who are founders with a single offspring (and
// no multiple matings) or who have no descendants
for (int i = fam->count - 1; i >= 0; i--)
{
if (discardable[i]) continue;
Person & p = *(persons[fam->path[i]]);
if (p.isFounder() || discardable[p.father->traverse])
{
if (mates[i] == -2 ||
offspring[i] > 1 ||
(mates[i] >= fam->founders &&
!discardable[persons[fam->path[mates[i]]]->father->traverse]) ||
haveData[p.traverse] ||
(mates[i] != -1 && haveData[mates[i]]))
continue;
trimming = true;
discardable[i] = true;
continue;
}
}
}
for (int i = fam->count - 1; i >= 0; i--)
if (discardable[i])
{
if (!quiet)
{
ShowTrimHeader(showHeader);
printf(" Removing person %s->%s: No data\n",
(const char *) fam->famid,
(const char *) persons[fam->path[i]]->pid);
}
delete persons[fam->path[i]];
}
else
newPersons[newCount++] = persons[fam->path[i]];
}
if (!showHeader)
printf("\n");
delete [] persons;
persons = newPersons;
count = newCount;
Sort();
}
Family
Definition PedigreeFamily.h:28
IntArray
Definition IntArray.h:24
Person
Definition PedigreePerson.h:32