http://cal.vet.upenn.edu/projects/genetic/inbreed/coeff/page1.htm
http://cal.vet.upenn.edu/projects/genetic/inbreed/coeff/page2.htm
Fx = the coefficient of inbreeding of individual x
n = number of generations from the sire (father) of x back to some ancestor common to the sire and dam
n' = number of generations from the dam (mother) of x back to some ancestor common to the sire and dam
Fa = coefficient of inbreeding of the common ancestor when that animal is itself inbred
∑ = summation of the separate contributions of each common ancestor (the summation here refers to all possible pathways by which an allele could be passed from the common ancestor to both parents)
Let's start with an example.
We have two common ancestors here: 1 and 2
First calculate the contribution from common ancestor 1.
n = the number of generations from the sire (in this case 3) to the common ancestor 1, that is 1
n' = the number of generations from the sire (in this case 3) to the common ancestor 1, that is also 1
(1/2)^(1+1+1)*(1+0) = 1/8.
Why Fa = 0? Because the inbreeding coefficient of common ancestor 1 is 0.
Contribution from common ancestor 2 is apparently also 1/8 due to symmetry.
Thus, the inbreeding coefficient of x Fx = 1/8 + 1/8 = 1/4 = 0.25
I decide to use Javascript to implement this equation.
A crucial task is to find all the common ancestors. In this case it is easy to spot 1 and 2. Things get trickier when father mates daughter. We may end up with some common ancestors which should be excluded from the calculation. Check this complex scenario.
In this case, by definition, F, C, E, D, G are all common ancestors of A and B, since B is A's daughter. In fact, F, C, E, D have no contribution towards inbreeding coefficient. However, G does make contribution because it can reach B through another path that doesn't contain A (GHB).
We break the task down into finding the ancestors of father and mother respectively, then find the common ones.
1. Not ony do we need to know the ancestor, but also we need to know the path. (In my program, I call it 'chain')
So for father A, we have 6 ancestors (including himself) but 7 paths
A
A C
A C F
A C E
A D
A D E
A D G
For mother B, we have 9 ancestors (everyone except x) but 11 paths
B
B A
B A C
B A C F
B A C E
B A D
B A D E
B A D G
B H
B H G
B H M
2. Compare each pair of paths between A and B. The last element of the path is a common ancestor candidate. Needless to say, the last element of the pair needs to match one another to start with. Then we check in addtion to the last element, is there any other common element? If there is, we need to ditch this common ancestor (e.g. E in ACE and BACE) too because it doesn't contribute towards inbreeding coefficient as discussed earlier. Finally we identify two common ancestors: A (A, AB) and G (ADG, BHG). The length of two paths combined minus 1 is actually the n+n'+1 part in the equation.
Technical note. For the first time, I try to use a Javascript unit test framework -- QUnit. It has similar syntax to JUnit, and I found it very pleasant to use.
Highlight of new stuff learnt.
- QUnit
- Javascript closure/module
- slice() method - used to copy array and get the last item of the array
| Scenario | Inbreeding coefficient |
|---|---|
| Full siblings | 0.25 |
| Father daughter | 0.25 |
| Half siblings | 0.125 |
| First cousins | 0.0625 |
| 3 generation full siblings | 0.5 |
Finally the code (https://github.com/sunmingtao/inbreed)
The Person.inbreedingCoefficient.complex test case corresponds to the 2nd complex scenario above.
inbreeding.js
'use strict';
function Person(name, father, mother){
var myFather, myMother;
doFather(father);
doMother(mother);
function doFather(father){
if (notNull(father)){
myFather = father;
return;
}
return myFather;
}
function doMother(mother){
if (notNull(mother)){
myMother = mother;
return;
}
return myMother;
}
function doHasFather(){
return doHasParent(myFather);
}
function doHasMother(){
return doHasParent(myMother);
}
function doHasParent(parent){
return notNull(parent);
}
function doInbreedingCoefficient(){
var fatherChains = doFatherChains();
logChains(fatherChains, "Father chains of "+name);
var motherChains = doMotherChains();
logChains(motherChains, "Mother chains of "+name);
var result = 0;
for (var i=0; i<fatherChains.length; i++){
for (var j=0; j<motherChains.length; j++){
var num = numberGenerationToCommonAncestor(fatherChains[i], motherChains[j]);
if (num) { //is a common ancestor
var commonAncestor = fatherChains[i].fathestAncestor();
console.log("Found common ancestor: "+commonAncestor.name+". Num is "+num);
result += Math.pow(0.5, num)*(1+commonAncestor.inbreedingCoefficient());
}
}
}
return result;
}
function doFatherChains(){
return doParentChains(myFather);
}
function doMotherChains(){
return doParentChains(myMother);
}
function doParentChains(parent){
var chains = [];
if (doHasParent(parent)){
chains.push(Chain([parent])); //parent itself
var ancestors = parent.chains(); //All parents' ancestors
for (var i = 0; i < ancestors.length; i++){
var chain = ancestors[i];
chain.prepend(parent);
chains.push(chain);
}
}
return chains;
}
function doChains(){
return doFatherChains().concat(doMotherChains());
}
var publicAPI = {
name: name,
father: doFather,
mother: doMother,
hasFather: doHasFather,
hasMother: doHasMother,
fatherChains: doFatherChains,
motherChains: doMotherChains,
chains: doChains,
inbreedingCoefficient: doInbreedingCoefficient
};
return publicAPI;
}
function Chain(persons){
var myChain = persons;
function getChain(){
return myChain;
}
function doPrepend(person){
myChain.unshift(person);
}
function getLength(){
return myChain.length;
}
function getFathestAncestor(){
return myChain.slice(-1)[0];
}
function doToString(){
var result = '';
for (var i=0; i<persons.length; i++){
result += persons[i].name + ' ';
}
return result;
}
return {
chain : getChain,
prepend: doPrepend,
length: getLength,
fathestAncestor: getFathestAncestor,
toString: doToString
};
}
//returns the number of nodes to go from start point of chain 1 (item 0) to common ancestor, then to the start point of chain 2
//if two pathways to the common ancestor share a portion of same route, return 0. chain 1 = [a, c, f, g], chain 2 = [b, d, f, g].
//In this case, f is deemed 'common ancestor', g is not (because fg is shared by both chains)
//e.g.1 chain 1 = [a, c, f], chain 2 = [b, f], then traverse acfb, return 4
//e.g.2 the start point itself is a common ancestor. chain 1 = [a], chain 2 = [b, a], then traverse ab, return 2
//e.g.3 if no common ancestor, return 0. chain 1 = [a, c], chain 2 = [b, d]
function numberGenerationToCommonAncestor(chain1, chain2){
if (chain1.fathestAncestor() != chain2.fathestAncestor() || hasCloserCommonAncestor(chain1, chain2)){
return 0;
}
return chain1.length() + chain2.length() - 1;
}
//Check if other closer common ancestor exists apart from the farthest common ancestor
function hasCloserCommonAncestor(chain1, chain2){
var persons1 = chain1.chain().slice(); //make a copy
var persons2 = chain2.chain().slice();
persons1.pop(); //chop off the last item
persons2.pop();
for (var i=0; i<persons1.length; i++){
for (var j=0; j<persons2.length; j++){
if (persons1[i] == persons2[j]){
return true;
}
}
}
return false;
}
function logChains(chains, message){
console.log(message);
for (var i=0; i<chains.length; i++) {
console.log(chains[i].toString());
}
}
function notNull(v){
return typeof v !== 'undefined' && v !== null;
}
tests.js
QUnit.test( "Person.init", function( assert ) {
var a = Person('a');
assert.equal( a.name, 'a');
var b = Person('b');
var c = Person('c', a , b);
assert.equal(c.mother().name, b.name);
var d = Person('d');
d.father(a);
d.mother(b);
assert.equal(d.mother().name, b.name);
});
QUnit.test( "Person.hasFatherMother", function( assert ) {
var a = Person('a');
assert.equal( a.hasFather(), false);
var b = Person('b');
var c = Person('c', a , b);
assert.equal( c.hasFather(), true);
var d = Person('d');
d.father(a);
d.mother(b);
assert.equal( d.hasMother(), true);
});
QUnit.test( "Person.chain", function( assert ) {
var a = Person('a');
assert.equal( a.fatherChains().length, 0);
var b = Person('b');
var c = Person('c', a , b);
assert.equal( c.fatherChains().length, 1);
assert.equal( c.fatherChains()[0].chain()[0].name, 'a');
});
QUnit.test( "Person.longer.chain", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c', a , b);
var d = Person('d');
var e = Person('e', c, d);
assert.equal( e.fatherChains().length, 3);
assert.equal( e.fatherChains()[0].chain().length, 1);
assert.equal( e.fatherChains()[0].chain()[0].name, 'c');
assert.equal( e.fatherChains()[1].chain().length, 2);
assert.equal( e.fatherChains()[1].chain()[0].name, 'c');
assert.equal( e.fatherChains()[1].fathestAncestor().name, 'a');
assert.equal( e.fatherChains()[2].chain().length, 2);
assert.equal( e.fatherChains()[2].chain()[0].name, 'c');
assert.equal( e.fatherChains()[2].fathestAncestor().name, 'b');
});
QUnit.test( "Chain.prepend", function( assert ) {
var a = Person('a');
var chain = Chain([a]);
assert.equal( chain.length(), 1);
var b = Person('b');
chain.prepend(b);
assert.equal( chain.length(), 2);
assert.equal( chain.chain()[0].name, 'b');
});
QUnit.test("numberGenerationToCommonAncestor", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c');
var d = Person('d');
var e = Person('e');
var chain1 = Chain([a]);
var chain2 = Chain([b]);
assert.equal(numberGenerationToCommonAncestor(chain1, chain2), 0);
var chain1 = Chain([c, b]);
var chain2 = Chain([b]);
assert.equal(numberGenerationToCommonAncestor(chain1, chain2), 2);
var chain1 = Chain([c, a]);
var chain2 = Chain([d, a]);
assert.equal(numberGenerationToCommonAncestor(chain1, chain2), 3);
var chain1 = Chain([c, a, b]);
var chain2 = Chain([d, a, b]);
assert.equal(numberGenerationToCommonAncestor(chain1, chain2), 0);
});
QUnit.test("Person.inbreedingCoefficient.notInbred", function( assert ) {
var a = Person('a');
var b = Person('b');
var x = Person('x', a, b);
assert.equal(x.inbreedingCoefficient(), 0);
});
QUnit.test("Person.inbreedingCoefficient.fullSibling", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c', a, b);
var d = Person('d', a, b);
var x = Person('x', c, d);
assert.equal(x.inbreedingCoefficient(), 0.25);
});
QUnit.test("Person.inbreedingCoefficient.halfSibling", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c');
var d = Person('d', a, b);
var e = Person('e', b, c);
var x = Person('x', d, e);
assert.equal(x.inbreedingCoefficient(), 0.125);
});
QUnit.test("Person.inbreedingCoefficient.fatherDaught", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c', a, b);
var x = Person('x', a, c);
assert.equal(x.inbreedingCoefficient(), 0.25);
});
QUnit.test("Person.inbreedingCoefficient.firstCousin", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c', a, b);
var d = Person('d', a, b);
var e = Person('e');
var f = Person('f');
var g = Person('g', c, e);
var h = Person('h', d, f);
var x = Person('x', g, h);
assert.equal(x.inbreedingCoefficient(), 0.0625);
});
QUnit.test("Person.inbreedingCoefficient.complex", function( assert ) {
var F = Person('F');
var E = Person('E');
var G = Person('G');
var M = Person('M');
var C = Person('C', F, E);
var D = Person('D', E, G);
var H = Person('H', G, M);
var A = Person('A', C, D);
var B = Person('B', A, H);
var X = Person('X', A, B);
assert.equal(X.inbreedingCoefficient(), 0.3125);
});
QUnit.test("Person.inbreedingCoefficient.deepInbred", function( assert ) {
var a = Person('a');
var b = Person('b');
var c = Person('c', a, b);
var d = Person('d', a, b);
var e = Person('e', c, d);
var f = Person('f', c, d);
var g = Person('g', e, f);
var h = Person('h', e, f);
var x = Person('x', g, h);
assert.equal(x.inbreedingCoefficient(), 0.5);
});
html
<!DOCTYPE html>
<html>
<head>
<meta charset="utf-8">
<title>QUnit basic example</title>
<link rel="stylesheet" href="https://code.jquery.com/qunit/qunit-2.0.1.css">
</head>
<body>
<div id="qunit"></div>
<div id="qunit-fixture"></div>
<script src="https://code.jquery.com/qunit/qunit-2.0.1.js"></script>
<script type="text/javascript" src="inbreeding.js"></script>
<script type="text/javascript" src="tests.js"></script>
</body>
</html>
