deps/kb_meme/meme-5.0.1/src/tomtom.c
/**************************************************************************
 * Estimate the source distribution from counts, using a 1-component
 * or 5-component Dirichlet prior.
 *
 * The 5-component prior only works for DNA.
 **************************************************************************/
Estimate the source distribution from counts, using a 1-component
or 5-component Dirichlet prior.

The 5-component prior only works for DNA.
Estimate the source distribution from counts using a 1 - component or 5 - component Dirichlet prior . The 5 - component prior only works for DNA .
static ARRAY_T* estimate_source
(int      num_dirichlets,
 ARRAY_T* counts)
{
  int dna_hypers[6][4] = {
    {1, 1, 1, 1}, 
    {5, 1, 1, 1}, 
    {1, 5, 1, 1},
    {1, 1, 5, 1},
    {1, 1, 1, 5},
    {2, 2, 2, 2}
  };
  int protein_hypers[20] 
    = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};
   or 5-component mixture?
  int alph_size = get_array_length(counts);
  if (num_dirichlets == 1) {
    if (alph_size == 4) {
      return(one_dirichlet_component(dna_hypers[0], counts));
    } else if (alph_size == 20) {
      return(one_dirichlet_component(protein_hypers, counts));
    } else {
      die("Invalid alphabet size (%d) for BLiC score.\n", alph_size);
    }
  }
  if (alph_size != 4) {
    die("Sorry, BLiC5 is only implemented for DNA motifs.\n");
  }
  int i_hyper;
  double loga=log(1.0/(float) num_dirichlets)+compute_log_likelihood(dna_hypers[1], counts);
  double logb = 0.0;
  double log_sum = 0.0;
  for (i_hyper = 2; i_hyper <= num_dirichlets; i_hyper++) {
    logb=log(1.0/(float) num_dirichlets)+compute_log_likelihood(dna_hypers[i_hyper], counts);
    log_sum = logsum(loga,logb);
    loga = log_sum; 
  }
  double log_denominator=log_sum;
  ARRAY_T* return_value = allocate_array(alph_size);
  for (i_hyper = 1; i_hyper <= num_dirichlets; i_hyper++) {
    double log_posterior = log(1.0 / (float) num_dirichlets) 
      + compute_log_likelihood(dna_hypers[i_hyper], counts) - log_denominator;
    ARRAY_T* one_component
      = one_dirichlet_component(dna_hypers[i_hyper], counts);
    scalar_mult(exp(log_posterior), one_component);
    sum_array(one_component, return_value);
    free_array(one_component);
  }
  return(return_value);
}
static ARRAY_T * estimate_source ( int num_dirichlets , ARRAY_T * counts ) { int dna_hypers [ 6 ] [ 4 ] = { { 1 , 1 , 1 , 1 } , { 5 , 1 , 1 , 1 } , { 1 , 5 , 1 , 1 } , { 1 , 1 , 5 , 1 } , { 1 , 1 , 1 , 5 } , { 2 , 2 , 2 , 2 } } ; int protein_hypers [ 20 ] = { 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 } ; int alph_size = get_array_length ( counts ) ; if ( num_dirichlets == 1 ) { if ( alph_size == 4 ) { return ( one_dirichlet_component ( dna_hypers [ 0 ] , counts ) ) ; } else if ( alph_size == 20 ) { return ( one_dirichlet_component ( protein_hypers , counts ) ) ; } else { die ( " \n " , alph_size ) ; } } if ( alph_size != 4 ) { die ( " \n " ) ; } int i_hyper ; double loga = log ( 1.0 / ( float ) num_dirichlets ) + compute_log_likelihood ( dna_hypers [ 1 ] , counts ) ; double logb = 0.0 ; double log_sum = 0.0 ; for ( i_hyper = 2 ; i_hyper <= num_dirichlets ; i_hyper ++ ) { logb = log ( 1.0 / ( float ) num_dirichlets ) + compute_log_likelihood ( dna_hypers [ i_hyper ] , counts ) ; log_sum = logsum ( loga , logb ) ; loga = log_sum ; } double log_denominator = log_sum ; ARRAY_T * return_value = allocate_array ( alph_size ) ; for ( i_hyper = 1 ; i_hyper <= num_dirichlets ; i_hyper ++ ) { double log_posterior = log ( 1.0 / ( float ) num_dirichlets ) + compute_log_likelihood ( dna_hypers [ i_hyper ] , counts ) - log_denominator ; ARRAY_T * one_component = one_dirichlet_component ( dna_hypers [ i_hyper ] , counts ) ; scalar_mult ( exp ( log_posterior ) , one_component ) ; sum_array ( one_component , return_value ) ; free_array ( one_component ) ; } return ( return_value ) ; }
Estimate the source distribution from counts, using a 1-component
 or 5-component Dirichlet prior.
Estimate the source distribution from counts using a 1 - component or 5 - component Dirichlet prior .
// Single component // Five components ... // Single component or 5-component mixture? // Die if trying to do proteins. // Compute the denominator. /* 
       Compute the posterior.
       \Pr(k|n) = \frac{\Pr(k) \Pr(n|k)}{\sum_j \Pr(j) \Pr(n|j)}
       */ // Compute the Dirichlet component. // Multiply them together and add to the return value.