tdg09.mli 4.6 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11
open Phylogenetics

module Evolution_model : sig
  type param = {
    stationary_distribution : Amino_acid.vector ;
    exchangeability_matrix : Amino_acid.matrix ;
    scale : float ;
  }
  val param_of_wag : Wag.t -> float -> param
  val rate_matrix : param -> Amino_acid.matrix
  val stationary_distribution : param -> Amino_acid.vector
12
  val transition_probability_matrix : param -> float -> Linear_algebra.Lacaml.mat
13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118
end

type likelihood_ratio_test = {
  full_log_likelihood : float ;
  reduced_log_likelihood : float ;
  _D_ : float ;
  df : float ;
  pvalue : float ;
}

module type S = sig
  type branch_info
  type leaf_info
  type site

  type simulation = (Amino_acid.t, Amino_acid.t, branch_info) Tree.t

  module Model1 : sig
    type param = float

    val maximum_log_likelihood :
      ?debug:bool ->
      exchangeability_matrix:Rate_matrix.Amino_acid.t ->
      stationary_distribution:Amino_acid.vector ->
      (_, leaf_info, branch_info) Tree.t ->
      site ->
      float * param

    val simulate_site :
      exchangeability_matrix:Amino_acid.matrix ->
      stationary_distribution:Amino_acid.vector ->
      (_, _, branch_info) Tree.t ->
      param:param ->
      simulation
  end

  module Model2 : sig
    type param = {
      scale : float ;
      stationary_distribution : Amino_acid.vector ;
    }

    val maximum_log_likelihood :
      ?debug:bool ->
      ?mode:[< `dense | `sparse > `sparse ] ->
      exchangeability_matrix:Rate_matrix.Amino_acid.t ->
      (_, leaf_info, branch_info) Tree.t ->
      site ->
      float * param

    val lrt :
      ?mode:[< `dense | `sparse > `sparse ] ->
      Wag.t ->
      (_, leaf_info, branch_info) Tree.t ->
      site ->
      Model1.param * param * likelihood_ratio_test
  end

  module Model3 : sig
    type param = {
      scale : float ;
      stationary_distribution0 : Amino_acid.vector ;
      stationary_distribution1 : Amino_acid.vector ;
    }

    val maximum_log_likelihood :
      ?debug:bool ->
      ?mode:[< `dense | `sparse > `sparse ] ->
      exchangeability_matrix:Rate_matrix.Amino_acid.t ->
      (_, leaf_info, branch_info) Tree.t ->
      site ->
      float * param

    val lrt :
      ?mode:[< `dense | `sparse > `sparse ] ->
      Wag.t ->
      (_, leaf_info, branch_info) Tree.t ->
      site ->
      Model2.param * param * likelihood_ratio_test

    val simulate_site :
      exchangeability_matrix:Rate_matrix.Amino_acid.t ->
      scale:float ->
      stationary_distribution0:Amino_acid.vector ->
      stationary_distribution1:Amino_acid.vector ->
      (_, leaf_info, branch_info) Tree.t ->
      simulation
  end
end

module type Leaf_info = sig
  type t
  type species
  val species : t -> species
  val condition : t -> [`Ancestral | `Convergent]
end

module type Branch_info = sig
  type t
  val length : t -> float
  val condition : t -> [`Ancestral | `Convergent]
end

module type Site = sig
  type t
  type species
Philippe Veber's avatar
Philippe Veber committed
119
  val get_aa : t -> species -> Amino_acid.t option
120 121 122 123 124 125 126
end

module Make(Branch_info : Branch_info)(Leaf_info : Leaf_info)(Site : Site with type species = Leaf_info.species) :
  S with type site := Site.t
     and type leaf_info := Leaf_info.t
     and type branch_info := Branch_info.t

127 128 129 130 131 132 133 134 135 136
module Pack : sig
  type leaf_info = int * Convergence_tree.condition

  module Leaf_info : Leaf_info with type t = leaf_info
                                and type species = int

  module Site : Site with type species = int
                      and type t = Amino_acid.t array

  include module type of Make(Convergence_tree.Branch_info)(Leaf_info)(Site)
137 138 139

  val simulate_profile : float -> Amino_acid.vector

140 141 142 143 144
  val pair_tree :
    branch_length1:float ->
    branch_length2:float ->
    npairs:int ->
    (unit, leaf_info, Convergence_tree.branch_info) Tree.t
145 146 147 148

  val convergence_tree :
    Convergence_tree.t ->
    (unit, leaf_info, Convergence_tree.branch_info) Tree.t
149 150 151 152 153
end

module Implementation_check : sig
  open Pack

154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183
  val likelihood_plot_demo : Wag.t -> unit

  val lrt_1_vs_2_null_simulation :
    ?seed:int ->
    ?mode:[< `dense | `sparse > `sparse ] ->
    ?nb_simulations:int ->
    Wag.t ->
    (simulation * Model1.param * Model2.param * likelihood_ratio_test) array

  val lrt_2_vs_3_null_simulation :
    ?seed:int ->
    ?mode:[< `dense | `sparse > `sparse ] ->
    ?alpha:float ->
    ?nb_simulations:int ->
    Wag.t ->
    (simulation * Model2.param * Model3.param * likelihood_ratio_test) array

  val render_pvalue_histogram :
    title:string ->
    (_ * _ * _ * likelihood_ratio_test) array ->
    string ->
    unit

  val render_stat_histogram :
    title:string ->
    df:float ->
    (_ * _ * _ * likelihood_ratio_test) array ->
    string ->
    unit
end