///|
/// Maximum parsimony tree construction demo (Biopython
/// ``Bio.Phylo.TreeConstruction`` equivalents): ParsimonyScorer (Fitch and
/// Sankoff), NNITreeSearcher hill climbing, ParsimonyTreeConstructor with
/// the default identity-UPGMA start, and Tree.root_at_midpoint.
///
/// Run: moon run examples/parsimony_tree_demo
fn main raise {
  println("=== Bio.Phylo.TreeConstruction Parsimony Demo ===")
  println("")

  // -------------------------------------------------------------------------
  // 1. Five-taxon nucleotide alignment
  // -------------------------------------------------------------------------
  let pairs = [
    ("Alpha", "AACGTGGCCACAT"),
    ("Beta", "AAGGTCGCCACAC"),
    ("Gamma", "CAGTTCGCCACAA"),
    ("Delta", "GAGATTTCCGCCT"),
    ("Epsilon", "GAGATCTCCGCCC"),
  ]
  let alignment = make_map(pairs)
  println("--- 1. Input alignment (5 taxa, 13 columns) ---")
  for p in pairs {
    println("  \{p.0}: \{p.1}")
  }
  println("")

  // Neighbor-joining start tree (unrooted), as built by Biopython.
  let calculator = @bio.DistanceCalculator::new(model="identity")
  let dm = calculator.get_distance(pairs)
  let names = pairs.map(fn(kv) { kv.0 })
  let nj = @bio.Tree::new(
    @bio.TreeConstructor::new(tree_method="nj").build_tree(dm, names),
    rooted=false,
  )
  println("--- 2. Unrooted NJ start tree ---")
  println("  " + nj.root.to_newick())
  println("")

  // Fitch parsimony: the scorer roots unrooted trees at their midpoint.
  let fitch = @bio.ParsimonyScorer::new()
  let midpoint = nj.root_at_midpoint()
  println("--- 3. Automatic midpoint rooting before scoring ---")
  println("  " + midpoint.root.to_newick())
  println("  Fitch parsimony score: \{fitch.get_score(nj, alignment)}")
  println("")

  // NNI neighbors and the hill-climbing search.
  println("--- 4. NNI hill climbing (Fitch) ---")
  let neighbors = @bio.nni_get_neighbors(midpoint)
  let neighbor_scores : Array[String] = []
  for n in neighbors {
    neighbor_scores.push(fitch.get_score(n, alignment).to_string())
  }
  println("  neighbor scores: [" + neighbor_scores.join(",") + "]")
  let result = @bio.NNITreeSearcher::new(fitch).search(nj, alignment)
  println("  final tree: " + result.root.to_newick())
  println("  final score: \{fitch.get_score(result, alignment)}")
  println("")

  // Default constructor: identity-UPGMA start, then NNI search.
  println("--- 5. ParsimonyTreeConstructor (default UPGMA start) ---")
  let ctor = @bio.ParsimonyTreeConstructor::new(
    @bio.NNITreeSearcher::new(@bio.ParsimonyScorer::new()),
  )
  let built = ctor.build_tree(alignment)
  println("  " + built.root.to_newick())
  println("  Fitch score: \{fitch.get_score(built, alignment)}")
  println("")

  // -------------------------------------------------------------------------
  // 2. Six-taxon Sankoff parsimony (transitions cost 1, transversions 5)
  // -------------------------------------------------------------------------
  println("--- 6. Sankoff NNI search with a custom cost matrix ---")
  let sankoff_alignment = make_map([
    ("S1", "AAAAA"),
    ("S2", "AAGAA"),
    ("S3", "GTATT"),
    ("S4", "GTAAT"),
    ("S5", "TTTTT"),
    ("S6", "TTTGT"),
  ])
  let start = @bio.Tree::new(
    node(
      [
        node(
          [
            node(
              [
                node([leaf("S1", 0.11), leaf("S3", 0.13)], 0.31),
                leaf("S4", 0.14),
              ],
              0.32,
            ),
            leaf("S5", 0.15),
          ],
          0.33,
        ),
        node([leaf("S6", 0.16), leaf("S2", 0.12)], 0.34),
      ],
      0.02,
    ),
    rooted=true,
  )
  let matrix = @bio.ParsimonyMatrix::new(["A", "C", "G", "T"], [
    [0.0, 5.0, 1.0, 5.0],
    [5.0, 0.0, 5.0, 1.0],
    [1.0, 5.0, 0.0, 5.0],
    [5.0, 1.0, 5.0, 0.0],
  ])
  let sankoff = @bio.ParsimonyScorer::new(matrix=Some(matrix))
  println("  start: " + start.root.to_newick())
  println("  Fitch score:   \{fitch.get_score(start, sankoff_alignment)}")
  println("  Sankoff score: \{sankoff.get_score(start, sankoff_alignment)}")

  let mut work = start
  let mut round = 1
  let mut done = false
  while !done {
    let current = sankoff.get_score(work, sankoff_alignment)
    let round_neighbors = @bio.nni_get_neighbors(work)
    let scores : Array[String] = []
    let mut best_score = current
    let mut accepted = -1
    for i in 0.. Map[String, String] {
  let m : Map[String, String] = Map([])
  for p in pairs {
    m[p.0] = p.1
  }
  m
}

///|
fn leaf(name : String, bl : Double) -> @bio.Clade {
  @bio.Clade::new(name=Some(name), branch_length=Some(bl))
}

///|
fn node(children : Array[@bio.Clade], bl : Double) -> @bio.Clade {
  @bio.Clade::new(branch_length=Some(bl), clades=children)
}