eendeego/.gitignore

## .gitignore

      
    Raw
  

              .gitignore
            
          
## README.md

      
    Raw
  

              README.md
            
          
    This is a quick hack demonstrating a Moran Model.

  
## index.html
<html>
<head>
  <title>Moran Model Simulator</title>
<style type="text/css">
body {
  background: #aaa;
}

svg {
  background: #ccc;
}

* {
  font-family: Helvetica, sans;
}

pre {
  font-family: Menlo, Monaco, monospace;
}

label {
  display: inline-block;
  width: 8em;
}

#data {
  display: -moz-box;
  display: -webkit-box;
  display: box;
  flex-direction: lr;
  width: 100%;
}

#data form {
/*  display: inline-block;*/
  min-width: 50%;
  max-width: 50%;
  width: 50%;
}

#data div {
/*  display: inline-block;*/
  min-width: 50%;
  max-width: 50%;
  width: 50%;
}

#profile-graph circle.node {
  stroke: #fff;
  stroke-width: 1.5px;
}

#profile-graph line.link {
  stroke: #999;
  stroke-opacity: .6;
}

#sid-chart {
  margin-top: 2px;
}

#sid-chart path {
  stroke: 88f;
  stroke-width: 2;
  fill: none;
}

</style>
</head>
<body>
  <div id="profile-graph"></div>
  <div id="sid-chart"></div>
  <div id="data">
    <form>
      <label for="seed">Random seed: </label><input type="text" id="seed"><br>
      <label for="population">Population: </label><input type="text" id="population"><br>
      <label for="pmutation">P<sub>mutation</sub>: </label><input type="text" id="pmutation"><br>
      <label for="precomb">P<sub>recomb</sub>: </label><input type="text" id="precomb"><br>
      <label for="iterations">Iterations: </label><input type="text" id="iterations">
      <input type="button" value="start!" onclick="javascript:init(false);">
      <input type="button" value="start+iterate!" onclick="javascript:init(true);">
      <input type="button" value="iterate!" onclick="javascript:iterate();"><br>
      <input type="button" value="force.start!" onclick="javascript:force.start();">
      <input type="button" value="redraw!" onclick="javascript:redraw();">
      <input type="button" value="resize!" onclick="javascript:resize();">
    </form>
    <div id="stats">
      <p>Iterations: <span id="totalIterations"></span></p>
      <p>SID<sub>A</sub> = <span id="sid"></span></p>
      <p>CI<sub>95%</sub> = <span id="ci95"></span></p>
      <p>var<sub>SID</sub> = <span id="varsid"></span></p>
    </div>
  </div>
  <pre id="profiles"></pre>
<script type="text/javascript" src="//cdnjs.cloudflare.com/ajax/libs/d3/3.4.12/d3.min.js"></script>
<script type="text/javascript" src="random.js"></script>
<script type="text/javascript" src="individual.js"></script>
<script type="text/javascript" src="population.js"></script>
<script type="text/javascript" src="profile-tracker.js"></script>
<script type="text/javascript">

// based on http://stackoverflow.com/questions/2090551/parse-query-string-in-javascript
function parseQuery() {
  return window.location.search.substring(1).split("&").reduce(function(map, param) {
    var pair = param.split("=");
    map[pair[0]] = unescape(pair[1]);
    return map;
    }, {});
}

var population;
var profileTracker;

var w = 960,
    h = 500, h2 = 100,
    fill = d3.scale.category20();

var vis = d3.select("#profile-graph")
  .append("svg:svg")
    .attr("width", w)
    .attr("height", h);

var g_scaling = vis.append("svg:g");
var g_links   = g_scaling.append("svg:g").attr('id', 'links');
var g_nodes   = g_scaling.append("svg:g").attr('id', 'nodes');

var graph, nodes, links;

var force, link, node;

var sid_x, sid_y = d3.scale.linear().domain([0, h2]).range([0, h2]);

var sid_chart = d3.select("#sid-chart")
  .append("svg:svg")
    .attr("width", w)
    .attr("height", h2);

var sid_g = sid_chart.append("svg:g")
  .attr("transform", "translate(0, " + h2 + ")");

var sid_line = d3.svg.line()
    .x(function(d,i) { return sid_x(i); })
    .y(function(d) { return -100 * sid_y(d); })

function init(iterate) {
  var random_seed     = parseInt(document.getElementById("seed").value, 10);
  var population_size = parseInt(document.getElementById("population").value, 10);

  if(population !== undefined) {
    g_links.selectAll("line.link").remove();
    g_nodes.selectAll("circle.node").remove();
  }

  var random = new Random(random_seed);
  population = new Population(population_size, random);
  profileTracker = new ProfileTracker(population, population_size, random);

  population.initialize();
  if(iterate) { runIterations(); }

  showStats();

  graph = profileTracker.graph();
  nodes = graph.nodes;
  links = graph.links;

  force = d3.layout.force()
      .charge(-120)
      .linkDistance(30)
      .nodes(nodes)
      .links(links)
      .size([w, h]);

  force.on("tick", function() {
      link.attr("x1", function(d) { if(d.source.x!=d.source.x) { d.source.x=w/2;d.source.y=h/2; };return d.source.x; })
          .attr("y1", function(d) { return d.source.y; })
          .attr("x2", function(d) { if(d.target.x!=d.target.x) { d.target.x=w/2;d.target.y=h/2; };return d.target.x; })
          .attr("y2", function(d) { return d.target.y; });

      node.attr("cx", function(d) { if(d.x!=d.x) { d.x=w/2;d.y=h/2; };return d.x; })
          .attr("cy", function(d) { return d.y; });
    });

  redraw();

  force.start();

  vis.style("opacity", 1e-6)
    .transition()
      .duration(1000)
      .style("opacity", 1);
}


function redraw() {
  removed = []
  g_links.selectAll("line.link")
    .data(links, function(d) { return d.id; })
    .exit()
      .call(function(d) {
        removed.push(d.toString());
      })
      .remove();

  link = g_links.selectAll("line.link")
      .data(links, function(d) { return d.id; })
    .enter().append("svg:line")
      .attr("class", "link")
      .style("stroke-width", "1") // placeholder
      .style("stroke", "#fff")    // placeholder
      .attr("x1", function(d) { return d.source.x; })
      .attr("y1", function(d) { return d.source.y; })
      .attr("x2", function(d) { return d.target.x; })
      .attr("y2", function(d) { return d.target.y; });

  // all current links
  g_links.selectAll("line.link")
    .style("stroke-width", function(d) {
      if(d.value !== 1) {
        return d.value * 3;
      } else {
        return d.value * 3;
      };
    })
    .style("stroke", function(d) {
      if(d.mutation) {
        return d.hasRecombination ? 'hsl(120,100%,20%)' : 'hsl(0,100%,50%)';
      } else {
        return 'hsl(120,50%,' + Math.min(100, 40 + d.value * 10) + '%)';
      }
    });

  removed = []
  g_nodes.selectAll("circle.node")
      .data(nodes, function(d) { return d.id; })
    .exit()
      .call(function(d) {
        removed.push(d.toString());
      })
      .remove();

  // all current nodes
  g_nodes.selectAll("circle.node")
    .attr("r", function(d) { return 2 * Math.sqrt(d.count); })
    .transition()
    .duration(1000);

  // node => new nodes
  node = g_nodes.selectAll("circle.node")
      .data(nodes, function(d) { return d.id; })
    .enter().append("svg:circle")
      .attr("class", "node")
      .attr("cx", function(d) { return d.x; })
      .attr("cy", function(d) { return d.y; })
      .attr("r", function(d) { return 2 * Math.sqrt(d.count); })
      .style("fill", function(d) { return fill(d.group); })
      .call(force.drag);

  node.append("svg:title");

  // all current nodes
  titles = [];
  g_nodes.selectAll("circle.node")
    .selectAll("title")
    .text(function(d) { var title = d.id + " (#" + d.count + "): " + d.key; titles.push(title); return title; });

  link = g_links.selectAll("line.link");
  node = g_nodes.selectAll("circle.node");

  var sid_history = profileTracker.sidHistory();
  sid_x = d3.scale.linear().domain([0, Math.max(99, sid_history.length-1)]).range([0, w]);
  sid_g.selectAll("path").remove();
  sid_g.append("svg:path").attr("d", sid_line(sid_history));
}

function runIterations() {
  var mutationProbability = parseFloat(document.getElementById("pmutation").value);
  var recombProbability   = parseFloat(document.getElementById("precomb").value);
  var iterations          = parseInt(document.getElementById("iterations").value, 10);
  for(var i=0; i<iterations; i++) {
    population.generationStep(mutationProbability, recombProbability);
  }
}

function iterate() {
  runIterations();

  showStats();

  redraw();
  force.start();
}

function showStats() {
  var stats = profileTracker.stats();
  document.getElementById("totalIterations").innerHTML = stats.totalIterations;
  document.getElementById("sid").innerHTML = stats.SID.toPrecision(5);
  document.getElementById("ci95").innerHTML = "[" + stats.CI95.map(function(e) { return e.toPrecision(5); }).join(', ') + "]";
  document.getElementById("varsid").innerHTML = stats.varSID.toPrecision(5);

  document.getElementById("profiles").innerHTML = profileTracker.toCountsString();
}

function resize() {
  if(nodes === undefined || nodes.length == 0) {
    return;
  }

  var minx = nodes[0].x, miny = nodes[0].y,
      maxx = nodes[0].x, maxy = nodes[0].y;

  for(var i=0; i<nodes.length; i++) {
    var d = nodes[i];
    minx = Math.min(d.x,minx); maxx = Math.max(d.x,maxx);
    miny = Math.min(d.y,miny); maxy = Math.max(d.y,maxy);
  }

  minx = w/2 - minx; maxx = maxx - w/2;
  miny = h/2 - miny; maxy = maxy - h/2;

  maxx = Math.max(minx, maxx);
  maxy = Math.max(miny, maxy);

  var sx = 0.5 * w / maxx;
  var sy = 0.5 * h / maxy;
  var s = 0.9 * Math.min(sx, sy);

  g_scaling.attr("transform",
                 "translate(" + w/2 + "," + h/2 + ")" + " " +
                 "scale(" + s + ")" + " " +
                 "translate(-" + w/2 + ",-" + h/2 + ")");
}

(function() {
  var params = parseQuery();

  var random_seed         = parseInt  (params["seed"      ] ||    42);
  var population_size     = parseInt  (params["population"] ||   100);
  var mutationProbability = parseFloat(params["pmutation" ] || 0.008);
  var recombProbability   = parseFloat(params["precomb"   ] || 0.001);
  var iterations          = parseInt  (params["iterations"] ||  1000);

  document.getElementById("seed"      ).value = random_seed;
  document.getElementById("population").value = population_size;
  document.getElementById("pmutation" ).value = mutationProbability;
  document.getElementById("precomb"   ).value = recombProbability;
  document.getElementById("iterations").value = iterations;

  if(params.init !== undefined) {
    init(false);
  } else if(params.go !== undefined) {
    init(true);
  }
})();

</script>
</body>
</html>

## individual.js
var Individual = (function() {
  var populationCount = 0;

  return function(generation, gene, origin) {
    var id = ++populationCount;
    var key = gene.toString();

    return {
      id: function() { return id; },
      key: function() { return key; },
      generation: function() { return generation; },
      gene: function() { return gene; },
      origin: function() { return origin; },
      toString: function() {
        return "{Individual: #" + id + ", gene=[" + key + "]}";
      }
    };
  };
})();

## population.js
var Population = function(populationSize, random) {
  var GENE_SIZE = 7;

  if(random === undefined) {
    random = new Random();
  }

  var geneCounts = new Array(GENE_SIZE);
  for(var i=0; i<GENE_SIZE; i++) { geneCounts[i]=1; }

  var population = new Array(populationSize);
  var totalIterations = 0;

  var listeners = { birth: [], death: [], generation: [] };

  function notifyBirth(individual) {
    listeners.birth.map(function(callback) {
      callback(individual);
    });
  }

  function notifyDeath(individual) {
    listeners.death.map(function(callback) {
      callback(individual);
    });
  }

  function notifyGeneration() {
    listeners.generation.map(function(callback) {
      callback(totalIterations);
    });
  }

  function initialize() {
    var initialAllelicProfile = new Array(7);
    for(var i=0; i<GENE_SIZE; i++) { initialAllelicProfile[i]=1; }

    for(var i=0; i<populationSize; i++) {
      var individual = new Individual(totalIterations, initialAllelicProfile);
      population[i] = individual;
      notifyBirth(individual);
    }
  }

  function mutatingClone(individual, mutationProbability, recombinationProbability) {
    var key = individual.key();
    var gene = individual.gene();
    var newGene = gene;
    var originIndividuals = {};
    var origin = {
      clonedIndividual: key,
      individuals: originIndividuals,
      hasRecombination: false,
    };

    function cloneGene() {
      if(newGene === gene) {
        newGene = gene.slice(0); // Clone
      }
    }

    originIndividuals[key] = true;
    for(var i=0; i<GENE_SIZE; i++) {
      var god = random.nextFloat();
      if(god < mutationProbability) {
        cloneGene();
        newGene[i] = ++geneCounts[i];
      } else if(god < mutationProbability + recombinationProbability) {
        cloneGene();
        var sourceIndividual = population[random.nextIntCapped(populationSize)];
        newGene[i] = sourceIndividual.gene()[i];
        originIndividuals[sourceIndividual.key()] = false;
        origin.hasRecombination = true;
      }
    }
    return new Individual(totalIterations, newGene, origin);
  }

  var generationStep = function(mutationProbability,
                                recombinationProbability) {
    var dead    = random.nextIntCapped(populationSize);
    var doubled = random.nextIntCapped(populationSize - 1);
    if(doubled >= dead) { doubled++; }

    var newIndividual =
      mutatingClone(population[doubled],
                    mutationProbability,
                    recombinationProbability);
    notifyBirth(newIndividual);
    population.push(newIndividual);
    notifyDeath(population.splice(dead, 1)[0]);

    totalIterations++;
    notifyGeneration();
  };

  var on = function(event, callback) {
    if(event !== 'birth' && event !== 'death' && event !== 'generation') {
      throw new Exception('Unknown event: "' + event + '"');
    }

    listeners[event].push(callback);
  };

  var toString = function() {
    var result = "";
    for(var i=0; i<populationSize; i++) {
      result += population[i].toString() + "\n";
    }
    return result;
  }

  return {
    initialize: initialize,
    size: function() { return populationSize; },
    totalIterations: function() { return totalIterations; },
    generationStep: generationStep,
    on: on,
    toString: toString
  };
};

// var population = new Population(1000);
//
// console.log(population.toString());
//
// for(var i=0; i<1000; i++) {
//   population.generationStep(0.2);
// //  console.log(population.toString());
// }
//
// console.log(population.toString());

## profile-tracker.js
var ProfileTracker = function(population,
                              sidHistoryCycle,
                              random) {
  if(random === undefined) {
    random = new Random();
  }

  var profileInfoMap = {};
  var totalProfileCount = 0;

  var liveSet = [];
  var links = [];
  var linkMap = {};

  var sidHistory = [0];

  var linkId = function(source, target) {
    if(source < target) {
      return source + "-" + target;
    } else {
      return target + "-" + source;
    }
  }

  population.on('birth', function(individual) {
    var key = individual.key();
    var origin = individual.origin();
    var profileInfo = profileInfoMap[key];
    if(profileInfoMap[key] === undefined) {
      // New profile
      totalProfileCount++;
      profileInfo = profileInfoMap[key] = {
          key: key,
          gene: individual.gene(),
          count: 0,
          id: totalProfileCount,
          origin: origin,
          toString: function() { return this.id + " (#" + this.count + "): " + this.key; }
        };
      if(origin !== undefined) {
        parentNode = profileInfoMap[origin.clonedIndividual];
        var newAngle = random.nextFloat() * Math.PI * 2;
        var newDist = 10 + random.nextFloat() * 30;
        profileInfo.x = (parentNode.x || 0) + newDist * Math.cos(newAngle);
        profileInfo.y = (parentNode.y || 0) + newDist * Math.sin(newAngle);
      }
    } else {
      // Already existing profile
      if(origin !== undefined) {
        if(profileInfo.origin === undefined) {
          profileInfo.origin = { individuals: {} };
        }
        for(var k in origin.individuals) {
          profileInfo.origin.individuals[k] = origin.individuals[k];
        }
      }
    }
    if(++profileInfo.count === 1) {
      liveSet.push(profileInfo);
      if(origin !== undefined) {
        for(var k in origin.individuals) {
          var id = linkId(profileInfo.key, k);
          var link = linkMap[id];

          if(link === undefined) {
            link = { source: profileInfoMap[k],
                     target: profileInfo,
                     hasRecombination: origin.hasRecombination,
                     mutation: origin.individuals[k],
                     id: id,
                     value: 1 };
            linkMap[id] = link;
            links.push(link);
          } else {
            console.log("Already known: " + link.id);
            link.value++;
            //alert("aaaaaaaaaaaaa")
          }
        }
      }
    }
  });

  population.on('death', function(individual) {
    var profileInfo = profileInfoMap[individual.key()];
    if(--profileInfo.count === 0) {
      // TODO Optimize this (if at all possible)
      for(var i=0; i<liveSet.length; i++) {
        if(liveSet[i].id === profileInfo.id) {
          liveSet.splice(i, 1);
          break;
        }
      }

      for(var i=links.length-1; i>=0; i--) {
        if(links[i].source === profileInfo || links[i].target === profileInfo) {
          delete linkMap[links[i].id];
          links.splice(i, 1);
        }
      }
    }
  });

  if(sidHistoryCycle !== undefined) {
    population.on('generation', function(totalIterations) {
      if(totalIterations % sidHistoryCycle == 0) {
        sidHistory.push(sid());
      }
    });
  }

  var graph = function() {
    return { nodes: liveSet, links: links };
  };

  var sid = function() {
    var S = liveSet.length;
    var N = population.size();
    var NN1 = N * (N-1);
    var SID, s = 0;
    for(var i=0; i<S; i++) {
      var ni = liveSet[i].count;
      s += ni * (ni - 1)
    }
    return 1.0 - s * 1.0 / NN1;
  }

  var stats = function() {
    // http://darwin.phyloviz.net/ComparingPartitions/index.php?link=Tut4

    var S = liveSet.length;
    var N = population.size();
    var NN1 = N * (N-1);
    var SID = sid();

    var s1 = 0, s2 = 0;

    for(var i=0; i<S; i++) {
      var pi = liveSet[i].count * 1.0 / N;
      var pi2 = pi * pi;

      s1 += pi2;
      s2 += pi2 * pi;
    }

    var varSID = (4.0*N*(N-1)*(N-2)*s2 + 2.0*N*(N-1)*s1 + 2.0*N*(N-1)*(2*N-3)*s1*s1) / (NN1 * NN1);
    var doubleRootVarSID = 2*Math.sqrt(varSID);
    var CI95 = [SID - doubleRootVarSID, SID + doubleRootVarSID];

    return { SID: SID, varSID: varSID, CI95: CI95, totalIterations: population.totalIterations() };
  };

  var toCountsString = function() {
    var counts = [];
    for (var seq in profileInfoMap) {
      var profileInfo = profileInfoMap[seq];
      counts.push(profileInfo);
    }

    var result = counts.sort(function(a, b) {
      return b.count - a.count;
    });

    result = counts.reduce(function(buffer, profile) {
      if(profile.count !== 0) {
        //return buffer + profile.id + " (#" + profile.count + "): " + profile.key + "\n";
        return buffer + profile + "\n";
      } else {
        return buffer;
      }
    }, "");

    return result;
  }

  return {
    sidHistory: function() { return sidHistory; },
    graph: graph,
    toCountsString: toCountsString,
    stats: stats
  };
};

## random.js
Random = (function() {
  return function(seed) {
    // http://en.wikipedia.org/wiki/Random_number_generation
    var w, z;

    if(seed === undefined) {
      w = Math.floor(Math.random() * 65536);
      z = Math.floor(Math.random() * 65536);
    } else {
      z = seed >> 16;
      w = seed & 0xffff;
    }

    var next = function() {
      z = 36969 * (z & 65535) + (z >> 16);
      w = 18000 * (w & 65535) + (w >> 16);
      return (z << 16) + w; // 32-bit result
    }

    var nextInt = function() {
      return next();
    }

    var nextIntCapped = function(max_plus_1) {
      return (next() & 0x7fffffff) % max_plus_1;
    }

    var nextFloat = function() {
      return (next() & 0x7fffffff) * 1.0 / 0x3fffffff;
    }

    return {
      nextInt: nextInt,
      nextIntCapped: nextIntCapped,
      nextFloat: nextFloat
    };
  };
})();

// var rr = new Random();
// for(var i=0; i<10; i++) {
//   console.log("i: " + i + " : " + rr.nextInt());
// }
//
// for(var i=0; i<10; i++) {
//   console.log("i: " + i + " : " + rr.nextIntCapped(20));
// }
//
// for(var i=0; i<10; i++) {
//   console.log("i: " + i + " : " + rr.nextFloat());
// }
	<html>
	<head>
	<title>Moran Model Simulator</title>
	<style type="text/css">
	body {
	background: #aaa;
	}

	svg {
	background: #ccc;
	}

	* {
	font-family: Helvetica, sans;
	}

	pre {
	font-family: Menlo, Monaco, monospace;
	}

	label {
	display: inline-block;
	width: 8em;
	}

	#data {
	display: -moz-box;
	display: -webkit-box;
	display: box;
	flex-direction: lr;
	width: 100%;
	}

	#data form {
	/* display: inline-block;*/
	min-width: 50%;
	max-width: 50%;
	width: 50%;
	}

	#data div {
	/* display: inline-block;*/
	min-width: 50%;
	max-width: 50%;
	width: 50%;
	}

	#profile-graph circle.node {
	stroke: #fff;
	stroke-width: 1.5px;
	}

	#profile-graph line.link {
	stroke: #999;
	stroke-opacity: .6;
	}

	#sid-chart {
	margin-top: 2px;
	}

	#sid-chart path {
	stroke: 88f;
	stroke-width: 2;
	fill: none;
	}

	</style>
	</head>
	<body>
	<div id="profile-graph"></div>
	<div id="sid-chart"></div>
	<div id="data">
	<form>
	<label for="seed">Random seed: </label><input type="text" id="seed"><br>
	<label for="population">Population: </label><input type="text" id="population"><br>
	<label for="pmutation">P<sub>mutation</sub>: </label><input type="text" id="pmutation"><br>
	<label for="precomb">P<sub>recomb</sub>: </label><input type="text" id="precomb"><br>
	<label for="iterations">Iterations: </label><input type="text" id="iterations">
	<input type="button" value="start!" onclick="javascript:init(false);">
	<input type="button" value="start+iterate!" onclick="javascript:init(true);">
	<input type="button" value="iterate!" onclick="javascript:iterate();"><br>
	<input type="button" value="force.start!" onclick="javascript:force.start();">
	<input type="button" value="redraw!" onclick="javascript:redraw();">
	<input type="button" value="resize!" onclick="javascript:resize();">
	</form>
	<div id="stats">
	<p>Iterations: <span id="totalIterations"></span></p>
	<p>SID<sub>A</sub> = <span id="sid"></span></p>
	<p>CI<sub>95%</sub> = <span id="ci95"></span></p>
	<p>var<sub>SID</sub> = <span id="varsid"></span></p>
	</div>
	</div>
	<pre id="profiles"></pre>
	<script type="text/javascript" src="//cdnjs.cloudflare.com/ajax/libs/d3/3.4.12/d3.min.js"></script>
	<script type="text/javascript" src="random.js"></script>
	<script type="text/javascript" src="individual.js"></script>
	<script type="text/javascript" src="population.js"></script>
	<script type="text/javascript" src="profile-tracker.js"></script>
	<script type="text/javascript">

	// based on http://stackoverflow.com/questions/2090551/parse-query-string-in-javascript
	function parseQuery() {
	return window.location.search.substring(1).split("&").reduce(function(map, param) {
	var pair = param.split("=");
	map[pair[0]] = unescape(pair[1]);
	return map;
	}, {});
	}

	var population;
	var profileTracker;

	var w = 960,
	h = 500, h2 = 100,
	fill = d3.scale.category20();

	var vis = d3.select("#profile-graph")
	.append("svg:svg")
	.attr("width", w)
	.attr("height", h);

	var g_scaling = vis.append("svg:g");
	var g_links = g_scaling.append("svg:g").attr('id', 'links');
	var g_nodes = g_scaling.append("svg:g").attr('id', 'nodes');

	var graph, nodes, links;

	var force, link, node;

	var sid_x, sid_y = d3.scale.linear().domain([0, h2]).range([0, h2]);

	var sid_chart = d3.select("#sid-chart")
	.append("svg:svg")
	.attr("width", w)
	.attr("height", h2);

	var sid_g = sid_chart.append("svg:g")
	.attr("transform", "translate(0, " + h2 + ")");

	var sid_line = d3.svg.line()
	.x(function(d,i) { return sid_x(i); })
	.y(function(d) { return -100 * sid_y(d); })

	function init(iterate) {
	var random_seed = parseInt(document.getElementById("seed").value, 10);
	var population_size = parseInt(document.getElementById("population").value, 10);

	if(population !== undefined) {
	g_links.selectAll("line.link").remove();
	g_nodes.selectAll("circle.node").remove();
	}

	var random = new Random(random_seed);
	population = new Population(population_size, random);
	profileTracker = new ProfileTracker(population, population_size, random);

	population.initialize();
	if(iterate) { runIterations(); }

	showStats();

	graph = profileTracker.graph();
	nodes = graph.nodes;
	links = graph.links;

	force = d3.layout.force()
	.charge(-120)
	.linkDistance(30)
	.nodes(nodes)
	.links(links)
	.size([w, h]);

	force.on("tick", function() {
	link.attr("x1", function(d) { if(d.source.x!=d.source.x) { d.source.x=w/2;d.source.y=h/2; };return d.source.x; })
	.attr("y1", function(d) { return d.source.y; })
	.attr("x2", function(d) { if(d.target.x!=d.target.x) { d.target.x=w/2;d.target.y=h/2; };return d.target.x; })
	.attr("y2", function(d) { return d.target.y; });

	node.attr("cx", function(d) { if(d.x!=d.x) { d.x=w/2;d.y=h/2; };return d.x; })
	.attr("cy", function(d) { return d.y; });
	});

	redraw();

	force.start();

	vis.style("opacity", 1e-6)
	.transition()
	.duration(1000)
	.style("opacity", 1);
	}



	function redraw() {
	removed = []
	g_links.selectAll("line.link")
	.data(links, function(d) { return d.id; })
	.exit()
	.call(function(d) {
	removed.push(d.toString());
	})
	.remove();

	link = g_links.selectAll("line.link")
	.data(links, function(d) { return d.id; })
	.enter().append("svg:line")
	.attr("class", "link")
	.style("stroke-width", "1") // placeholder
	.style("stroke", "#fff") // placeholder
	.attr("x1", function(d) { return d.source.x; })
	.attr("y1", function(d) { return d.source.y; })
	.attr("x2", function(d) { return d.target.x; })
	.attr("y2", function(d) { return d.target.y; });

	// all current links
	g_links.selectAll("line.link")
	.style("stroke-width", function(d) {
	if(d.value !== 1) {
	return d.value * 3;
	} else {
	return d.value * 3;
	};
	})
	.style("stroke", function(d) {
	if(d.mutation) {
	return d.hasRecombination ? 'hsl(120,100%,20%)' : 'hsl(0,100%,50%)';
	} else {
	return 'hsl(120,50%,' + Math.min(100, 40 + d.value * 10) + '%)';
	}
	});

	removed = []
	g_nodes.selectAll("circle.node")
	.data(nodes, function(d) { return d.id; })
	.exit()
	.call(function(d) {
	removed.push(d.toString());
	})
	.remove();

	// all current nodes
	g_nodes.selectAll("circle.node")
	.attr("r", function(d) { return 2 * Math.sqrt(d.count); })
	.transition()
	.duration(1000);

	// node => new nodes
	node = g_nodes.selectAll("circle.node")
	.data(nodes, function(d) { return d.id; })
	.enter().append("svg:circle")
	.attr("class", "node")
	.attr("cx", function(d) { return d.x; })
	.attr("cy", function(d) { return d.y; })
	.attr("r", function(d) { return 2 * Math.sqrt(d.count); })
	.style("fill", function(d) { return fill(d.group); })
	.call(force.drag);

	node.append("svg:title");

	// all current nodes
	titles = [];
	g_nodes.selectAll("circle.node")
	.selectAll("title")
	.text(function(d) { var title = d.id + " (#" + d.count + "): " + d.key; titles.push(title); return title; });

	link = g_links.selectAll("line.link");
	node = g_nodes.selectAll("circle.node");

	var sid_history = profileTracker.sidHistory();
	sid_x = d3.scale.linear().domain([0, Math.max(99, sid_history.length-1)]).range([0, w]);
	sid_g.selectAll("path").remove();
	sid_g.append("svg:path").attr("d", sid_line(sid_history));
	}

	function runIterations() {
	var mutationProbability = parseFloat(document.getElementById("pmutation").value);
	var recombProbability = parseFloat(document.getElementById("precomb").value);
	var iterations = parseInt(document.getElementById("iterations").value, 10);
	for(var i=0; i<iterations; i++) {
	population.generationStep(mutationProbability, recombProbability);
	}
	}

	function iterate() {
	runIterations();

	showStats();

	redraw();
	force.start();
	}

	function showStats() {
	var stats = profileTracker.stats();
	document.getElementById("totalIterations").innerHTML = stats.totalIterations;
	document.getElementById("sid").innerHTML = stats.SID.toPrecision(5);
	document.getElementById("ci95").innerHTML = "[" + stats.CI95.map(function(e) { return e.toPrecision(5); }).join(', ') + "]";
	document.getElementById("varsid").innerHTML = stats.varSID.toPrecision(5);

	document.getElementById("profiles").innerHTML = profileTracker.toCountsString();
	}

	function resize() {
	if(nodes === undefined \|\| nodes.length == 0) {
	return;
	}

	var minx = nodes[0].x, miny = nodes[0].y,
	maxx = nodes[0].x, maxy = nodes[0].y;

	for(var i=0; i<nodes.length; i++) {
	var d = nodes[i];
	minx = Math.min(d.x,minx); maxx = Math.max(d.x,maxx);
	miny = Math.min(d.y,miny); maxy = Math.max(d.y,maxy);
	}

	minx = w/2 - minx; maxx = maxx - w/2;
	miny = h/2 - miny; maxy = maxy - h/2;

	maxx = Math.max(minx, maxx);
	maxy = Math.max(miny, maxy);

	var sx = 0.5 * w / maxx;
	var sy = 0.5 * h / maxy;
	var s = 0.9 * Math.min(sx, sy);

	g_scaling.attr("transform",
	"translate(" + w/2 + "," + h/2 + ")" + " " +
	"scale(" + s + ")" + " " +
	"translate(-" + w/2 + ",-" + h/2 + ")");
	}

	(function() {
	var params = parseQuery();

	var random_seed = parseInt (params["seed" ] \|\| 42);
	var population_size = parseInt (params["population"] \|\| 100);
	var mutationProbability = parseFloat(params["pmutation" ] \|\| 0.008);
	var recombProbability = parseFloat(params["precomb" ] \|\| 0.001);
	var iterations = parseInt (params["iterations"] \|\| 1000);

	document.getElementById("seed" ).value = random_seed;
	document.getElementById("population").value = population_size;
	document.getElementById("pmutation" ).value = mutationProbability;
	document.getElementById("precomb" ).value = recombProbability;
	document.getElementById("iterations").value = iterations;

	if(params.init !== undefined) {
	init(false);
	} else if(params.go !== undefined) {
	init(true);
	}
	})();

	</script>
	</body>
	</html>
	var Individual = (function() {
	var populationCount = 0;

	return function(generation, gene, origin) {
	var id = ++populationCount;
	var key = gene.toString();

	return {
	id: function() { return id; },
	key: function() { return key; },
	generation: function() { return generation; },
	gene: function() { return gene; },
	origin: function() { return origin; },
	toString: function() {
	return "{Individual: #" + id + ", gene=[" + key + "]}";
	}
	};
	};
	})();
	var Population = function(populationSize, random) {
	var GENE_SIZE = 7;

	if(random === undefined) {
	random = new Random();
	}

	var geneCounts = new Array(GENE_SIZE);
	for(var i=0; i<GENE_SIZE; i++) { geneCounts[i]=1; }

	var population = new Array(populationSize);
	var totalIterations = 0;

	var listeners = { birth: [], death: [], generation: [] };

	function notifyBirth(individual) {
	listeners.birth.map(function(callback) {
	callback(individual);
	});
	}

	function notifyDeath(individual) {
	listeners.death.map(function(callback) {
	callback(individual);
	});
	}

	function notifyGeneration() {
	listeners.generation.map(function(callback) {
	callback(totalIterations);
	});
	}

	function initialize() {
	var initialAllelicProfile = new Array(7);
	for(var i=0; i<GENE_SIZE; i++) { initialAllelicProfile[i]=1; }

	for(var i=0; i<populationSize; i++) {
	var individual = new Individual(totalIterations, initialAllelicProfile);
	population[i] = individual;
	notifyBirth(individual);
	}
	}

	function mutatingClone(individual, mutationProbability, recombinationProbability) {
	var key = individual.key();
	var gene = individual.gene();
	var newGene = gene;
	var originIndividuals = {};
	var origin = {
	clonedIndividual: key,
	individuals: originIndividuals,
	hasRecombination: false,
	};

	function cloneGene() {
	if(newGene === gene) {
	newGene = gene.slice(0); // Clone
	}
	}

	originIndividuals[key] = true;
	for(var i=0; i<GENE_SIZE; i++) {
	var god = random.nextFloat();
	if(god < mutationProbability) {
	cloneGene();
	newGene[i] = ++geneCounts[i];
	} else if(god < mutationProbability + recombinationProbability) {
	cloneGene();
	var sourceIndividual = population[random.nextIntCapped(populationSize)];
	newGene[i] = sourceIndividual.gene()[i];
	originIndividuals[sourceIndividual.key()] = false;
	origin.hasRecombination = true;
	}
	}
	return new Individual(totalIterations, newGene, origin);
	}

	var generationStep = function(mutationProbability,
	recombinationProbability) {
	var dead = random.nextIntCapped(populationSize);
	var doubled = random.nextIntCapped(populationSize - 1);
	if(doubled >= dead) { doubled++; }

	var newIndividual =
	mutatingClone(population[doubled],
	mutationProbability,
	recombinationProbability);
	notifyBirth(newIndividual);
	population.push(newIndividual);
	notifyDeath(population.splice(dead, 1)[0]);

	totalIterations++;
	notifyGeneration();
	};

	var on = function(event, callback) {
	if(event !== 'birth' && event !== 'death' && event !== 'generation') {
	throw new Exception('Unknown event: "' + event + '"');
	}

	listeners[event].push(callback);
	};

	var toString = function() {
	var result = "";
	for(var i=0; i<populationSize; i++) {
	result += population[i].toString() + "\n";
	}
	return result;
	}

	return {
	initialize: initialize,
	size: function() { return populationSize; },
	totalIterations: function() { return totalIterations; },
	generationStep: generationStep,
	on: on,
	toString: toString
	};
	};

	// var population = new Population(1000);
	//
	// console.log(population.toString());
	//
	// for(var i=0; i<1000; i++) {
	// population.generationStep(0.2);
	// // console.log(population.toString());
	// }
	//
	// console.log(population.toString());
	var ProfileTracker = function(population,
	sidHistoryCycle,
	random) {
	if(random === undefined) {
	random = new Random();
	}

	var profileInfoMap = {};
	var totalProfileCount = 0;

	var liveSet = [];
	var links = [];
	var linkMap = {};

	var sidHistory = [0];

	var linkId = function(source, target) {
	if(source < target) {
	return source + "-" + target;
	} else {
	return target + "-" + source;
	}
	}

	population.on('birth', function(individual) {
	var key = individual.key();
	var origin = individual.origin();
	var profileInfo = profileInfoMap[key];
	if(profileInfoMap[key] === undefined) {
	// New profile
	totalProfileCount++;
	profileInfo = profileInfoMap[key] = {
	key: key,
	gene: individual.gene(),
	count: 0,
	id: totalProfileCount,
	origin: origin,
	toString: function() { return this.id + " (#" + this.count + "): " + this.key; }
	};
	if(origin !== undefined) {
	parentNode = profileInfoMap[origin.clonedIndividual];
	var newAngle = random.nextFloat() * Math.PI * 2;
	var newDist = 10 + random.nextFloat() * 30;
	profileInfo.x = (parentNode.x \|\| 0) + newDist * Math.cos(newAngle);
	profileInfo.y = (parentNode.y \|\| 0) + newDist * Math.sin(newAngle);
	}
	} else {
	// Already existing profile
	if(origin !== undefined) {
	if(profileInfo.origin === undefined) {
	profileInfo.origin = { individuals: {} };
	}
	for(var k in origin.individuals) {
	profileInfo.origin.individuals[k] = origin.individuals[k];
	}
	}
	}
	if(++profileInfo.count === 1) {
	liveSet.push(profileInfo);
	if(origin !== undefined) {
	for(var k in origin.individuals) {
	var id = linkId(profileInfo.key, k);
	var link = linkMap[id];

	if(link === undefined) {
	link = { source: profileInfoMap[k],
	target: profileInfo,
	hasRecombination: origin.hasRecombination,
	mutation: origin.individuals[k],
	id: id,
	value: 1 };
	linkMap[id] = link;
	links.push(link);
	} else {
	console.log("Already known: " + link.id);
	link.value++;
	//alert("aaaaaaaaaaaaa")
	}
	}
	}
	}
	});

	population.on('death', function(individual) {
	var profileInfo = profileInfoMap[individual.key()];
	if(--profileInfo.count === 0) {
	// TODO Optimize this (if at all possible)
	for(var i=0; i<liveSet.length; i++) {
	if(liveSet[i].id === profileInfo.id) {
	liveSet.splice(i, 1);
	break;
	}
	}

	for(var i=links.length-1; i>=0; i--) {
	if(links[i].source === profileInfo \|\| links[i].target === profileInfo) {
	delete linkMap[links[i].id];
	links.splice(i, 1);
	}
	}
	}
	});

	if(sidHistoryCycle !== undefined) {
	population.on('generation', function(totalIterations) {
	if(totalIterations % sidHistoryCycle == 0) {
	sidHistory.push(sid());
	}
	});
	}

	var graph = function() {
	return { nodes: liveSet, links: links };
	};

	var sid = function() {
	var S = liveSet.length;
	var N = population.size();
	var NN1 = N * (N-1);
	var SID, s = 0;
	for(var i=0; i<S; i++) {
	var ni = liveSet[i].count;
	s += ni * (ni - 1)
	}
	return 1.0 - s * 1.0 / NN1;
	}

	var stats = function() {
	// http://darwin.phyloviz.net/ComparingPartitions/index.php?link=Tut4

	var S = liveSet.length;
	var N = population.size();
	var NN1 = N * (N-1);
	var SID = sid();

	var s1 = 0, s2 = 0;

	for(var i=0; i<S; i++) {
	var pi = liveSet[i].count * 1.0 / N;
	var pi2 = pi * pi;

	s1 += pi2;
	s2 += pi2 * pi;
	}

	var varSID = (4.0N(N-1)(N-2)s2 + 2.0N(N-1)s1 + 2.0N(N-1)(2N-3)s1s1) / (NN1 NN1);
	var doubleRootVarSID = 2*Math.sqrt(varSID);
	var CI95 = [SID - doubleRootVarSID, SID + doubleRootVarSID];

	return { SID: SID, varSID: varSID, CI95: CI95, totalIterations: population.totalIterations() };
	};

	var toCountsString = function() {
	var counts = [];
	for (var seq in profileInfoMap) {
	var profileInfo = profileInfoMap[seq];
	counts.push(profileInfo);
	}

	var result = counts.sort(function(a, b) {
	return b.count - a.count;
	});

	result = counts.reduce(function(buffer, profile) {
	if(profile.count !== 0) {
	//return buffer + profile.id + " (#" + profile.count + "): " + profile.key + "\n";
	return buffer + profile + "\n";
	} else {
	return buffer;
	}
	}, "");

	return result;
	}

	return {
	sidHistory: function() { return sidHistory; },
	graph: graph,
	toCountsString: toCountsString,
	stats: stats
	};
	};
	Random = (function() {
	return function(seed) {
	// http://en.wikipedia.org/wiki/Random_number_generation
	var w, z;

	if(seed === undefined) {
	w = Math.floor(Math.random() * 65536);
	z = Math.floor(Math.random() * 65536);
	} else {
	z = seed >> 16;
	w = seed & 0xffff;
	}

	var next = function() {
	z = 36969 * (z & 65535) + (z >> 16);
	w = 18000 * (w & 65535) + (w >> 16);
	return (z << 16) + w; // 32-bit result
	}

	var nextInt = function() {
	return next();
	}

	var nextIntCapped = function(max_plus_1) {
	return (next() & 0x7fffffff) % max_plus_1;
	}

	var nextFloat = function() {
	return (next() & 0x7fffffff) * 1.0 / 0x3fffffff;
	}

	return {
	nextInt: nextInt,
	nextIntCapped: nextIntCapped,
	nextFloat: nextFloat
	};
	};
	})();

	// var rr = new Random();
	// for(var i=0; i<10; i++) {
	// console.log("i: " + i + " : " + rr.nextInt());
	// }
	//
	// for(var i=0; i<10; i++) {
	// console.log("i: " + i + " : " + rr.nextIntCapped(20));
	// }
	//
	// for(var i=0; i<10; i++) {
	// console.log("i: " + i + " : " + rr.nextFloat());
	// }