deps/kb_meme/meme-5.0.1/src/mhmm.c
/***********************************************************************
 * Read the motifs, filter them, create an order and spacing if
 * possible and create the transition frequency matrix and the
 * spacer average matrix.
 ***********************************************************************/
Read the motifs, filter them, create an order and spacing if
possible and create the transition frequency matrix and the
spacer average matrix.
Read the motifs filter them create an order and spacing if possible and create the transition frequency matrix and the spacer average matrix .
void load_filter_process_motifs_for_hmm
(
  int request_n, 
  RBTREE_T *request_nums, 
  RBTREE_T *request_ids, 
  double request_ev, 
  double request_complexity, 
  const char *request_order, 
  double hit_pv, 
  const char *motif_filename, 
  const char *bg_source, 
  double motif_pseudo, 
  double trans_pseudo,  
  double spacer_pseudo, 
  bool keep_unused, 
  ALPH_T **alphabet, 
  RBTREE_T **motifs, 
  ARRAY_T **background, 
  ORDER_T **order_spacing, 
  MATRIX_T **transp_freq, 
  MATRIX_T **spacer_ave 
) {
  int (*filter)(void *context, int motif_num, MOTIF_T *motif);
  void *filter_context;
  char *filter_name;
  REQMOT_T requested;
  RBTREE_T *ids2nums; 
  RBNODE_T *node;
  MREAD_T *mread;
  MOTIF_T *motif;
  int motif_i, rc_motif_i, strands, i;
  bool created, has_reverse_strand;
  ARRAYLST_T *motif_occurrences;
  filter = NULL;
  filter_context = NULL;
  filter_name = "<no filter>";
  if (request_n != 0) {
    filter = select_top_n;
    filter_context = &request_n;
    filter_name = "-nmotifs";
  }
  if (rbtree_size(request_nums) != 0 || rbtree_size(request_ids) != 0) {
    if (filter != NULL) {
      die("Can't use -motif or -motif-id with %s.\n", filter_name);
    }
    requested.nums = request_nums;
    requested.ids = request_ids;
    filter = select_requested;
    filter_context = &requested;
    filter_name = "-motif or -motif-id";
  }
  if (request_ev != 0.0) {
    if (filter != NULL) {
      die("Can't use -ethresh with %s.\n", filter_name);
    }
    filter = select_better_e_value;
    filter_context = &request_ev;
  }
  if (request_complexity != 0.0) {
    if (filter != NULL) {
      die("Can't use -lowcomp with %s.\n", filter_name);
    }
    filter = select_more_complex;
    filter_context = &request_complexity;
  }
  if (request_order != NULL && filter != NULL) {
    die("Can't use -order with %s.\n", filter_name);
  }
  if (filter == NULL) filter = select_all;
  mread = mread_create(motif_filename, OPEN_MFILE | CALC_AMBIGS | SCANNED_SITES, true);
  mread_set_bg_source(mread, bg_source, NULL);
  mread_set_pseudocount(mread, motif_pseudo);
  has_reverse_strand = (mread_get_strands(mread) == 2);
  *motifs = rbtree_create(motif_num_cmp, rbtree_intcpy, free, NULL, destroy_motif);
  ids2nums = rbtree_create(rbtree_strcmp, rbtree_strcpy, free, rbtree_intcpy, free);
  for (motif_i = 1, rc_motif_i = -1; mread_has_motif(mread); ++motif_i, --rc_motif_i) {
    motif = mread_next_motif(mread);
    strands = filter(filter_context, motif_i, motif);
    if (strands & MOTIF_PLUS) rbtree_make(*motifs, &motif_i, motif);
    if (strands & MOTIF_MINUS) rbtree_make(*motifs, &rc_motif_i, dup_rc_motif(motif));
    if (strands) rbtree_make(ids2nums, get_motif_id(motif), &motif_i);
    else destroy_motif(motif);
  }
  *background = mread_get_background(mread);
  *alphabet = alph_hold(mread_get_alphabet(mread));
  motif_occurrences = mread_get_motif_occurrences(mread);
  mread_destroy(mread);
  *order_spacing = create_order_from_ids(request_order, ids2nums);
  rbtree_destroy(ids2nums);
  if (*order_spacing) {
    for (node = rbtree_first(*motifs); node != NULL; node = rbtree_next(node)) {
      motif_i = *((int*)rbtree_key(node));
      if (!order_contains(motif_i, *order_spacing)) {
        rbtree_delete(*motifs, node, NULL, NULL);
      }
    }
  }
  if (motif_occurrences != NULL && arraylst_size(motif_occurrences) > 0) {
    compute_transitions_and_spacers(
        *motifs,
        motif_occurrences,
        has_reverse_strand,
        hit_pv,
        order_spacing,
        transp_freq,
        spacer_ave);
  } else {
    compute_naive_transitions_and_spacers(
        rbtree_size(*motifs), 
        transp_freq, 
        spacer_ave);
  }
  normalize_spacer_counts(
      trans_pseudo, 
      spacer_pseudo,
      keep_unused,
      *transp_freq, 
      *spacer_ave);
  throw_out_unused_motifs(*motifs, *transp_freq, *spacer_ave);
  destroy_motif_occurrences(motif_occurrences);
}
void load_filter_process_motifs_for_hmm ( int request_n , RBTREE_T * request_nums , RBTREE_T * request_ids , double request_ev , double request_complexity , const char * request_order , double hit_pv , const char * motif_filename , const char * bg_source , double motif_pseudo , double trans_pseudo , double spacer_pseudo , bool keep_unused , ALPH_T * * alphabet , RBTREE_T * * motifs , ARRAY_T * * background , ORDER_T * * order_spacing , MATRIX_T * * transp_freq , MATRIX_T * * spacer_ave ) { int ( * filter ) ( void * context , int motif_num , MOTIF_T * motif ) ; void * filter_context ; char * filter_name ; REQMOT_T requested ; RBTREE_T * ids2nums ; RBNODE_T * node ; MREAD_T * mread ; MOTIF_T * motif ; int motif_i , rc_motif_i , strands , i ; bool created , has_reverse_strand ; ARRAYLST_T * motif_occurrences ; filter = NULL ; filter_context = NULL ; filter_name = " " ; if ( request_n != 0 ) { filter = select_top_n ; filter_context = & request_n ; filter_name = " " ; } if ( rbtree_size ( request_nums ) != 0 || rbtree_size ( request_ids ) != 0 ) { if ( filter != NULL ) { die ( " \n " , filter_name ) ; } requested . nums = request_nums ; requested . ids = request_ids ; filter = select_requested ; filter_context = & requested ; filter_name = " " ; } if ( request_ev != 0.0 ) { if ( filter != NULL ) { die ( " \n " , filter_name ) ; } filter = select_better_e_value ; filter_context = & request_ev ; } if ( request_complexity != 0.0 ) { if ( filter != NULL ) { die ( " \n " , filter_name ) ; } filter = select_more_complex ; filter_context = & request_complexity ; } if ( request_order != NULL && filter != NULL ) { die ( " \n " , filter_name ) ; } if ( filter == NULL ) filter = select_all ; mread = mread_create ( motif_filename , OPEN_MFILE | CALC_AMBIGS | SCANNED_SITES , true ) ; mread_set_bg_source ( mread , bg_source , NULL ) ; mread_set_pseudocount ( mread , motif_pseudo ) ; has_reverse_strand = ( mread_get_strands ( mread ) == 2 ) ; * motifs = rbtree_create ( motif_num_cmp , rbtree_intcpy , free , NULL , destroy_motif ) ; ids2nums = rbtree_create ( rbtree_strcmp , rbtree_strcpy , free , rbtree_intcpy , free ) ; for ( motif_i = 1 , rc_motif_i = -1 ; mread_has_motif ( mread ) ; ++ motif_i , -- rc_motif_i ) { motif = mread_next_motif ( mread ) ; strands = filter ( filter_context , motif_i , motif ) ; if ( strands & MOTIF_PLUS ) rbtree_make ( * motifs , & motif_i , motif ) ; if ( strands & MOTIF_MINUS ) rbtree_make ( * motifs , & rc_motif_i , dup_rc_motif ( motif ) ) ; if ( strands ) rbtree_make ( ids2nums , get_motif_id ( motif ) , & motif_i ) ; else destroy_motif ( motif ) ; } * background = mread_get_background ( mread ) ; * alphabet = alph_hold ( mread_get_alphabet ( mread ) ) ; motif_occurrences = mread_get_motif_occurrences ( mread ) ; mread_destroy ( mread ) ; * order_spacing = create_order_from_ids ( request_order , ids2nums ) ; rbtree_destroy ( ids2nums ) ; if ( * order_spacing ) { for ( node = rbtree_first ( * motifs ) ; node != NULL ; node = rbtree_next ( node ) ) { motif_i = * ( ( int * ) rbtree_key ( node ) ) ; if ( ! order_contains ( motif_i , * order_spacing ) ) { rbtree_delete ( * motifs , node , NULL , NULL ) ; } } } if ( motif_occurrences != NULL && arraylst_size ( motif_occurrences ) > 0 ) { compute_transitions_and_spacers ( * motifs , motif_occurrences , has_reverse_strand , hit_pv , order_spacing , transp_freq , spacer_ave ) ; } else { compute_naive_transitions_and_spacers ( rbtree_size ( * motifs ) , transp_freq , spacer_ave ) ; } normalize_spacer_counts ( trans_pseudo , spacer_pseudo , keep_unused , * transp_freq , * spacer_ave ) ; throw_out_unused_motifs ( * motifs , * transp_freq , * spacer_ave ) ; destroy_motif_occurrences ( motif_occurrences ) ; }
Read the motifs, filter them, create an order and spacing if
 possible and create the transition frequency matrix and the
 spacer average matrix.
Read the motifs filter them create an order and spacing if possible and create the transition frequency matrix and the spacer average matrix .
// IN select the top n positive strand motifs // IN select the requested motif numbers // IN select the requested motif ids // IN select motifs with a <= evalue // IN select motifs with a >= complexity // IN specifiy motifs and their ordering and spacing // IN specify the p value for keeping motif occurences // IN file containing motifs to load // IN source of background // IN pseudocount applied to motifs // IN pseudocount applied to transitions // IN pseudocount applied to spacers // IN keep the unused motifs // OUT alphabet of the motifs // OUT motifs with file position as key // OUT background from motif file // OUT order and spacing (uses motif key as identifiers) // OUT transition frequencies // OUT spacer averages // mapping of motif IDs to numbers /**********************************************
   * SELECT MOTIF FILTER
   **********************************************/ // filter top n // filter requested // filter e threshold // filter complexity threshold // filter order (but not when reading motifs) /**********************************************
   * READING THE MOTIFS
   **********************************************/ // Parse the order string.  // Remove unused motifs // Turn the raw motifs and motif occurrences into the elements of the model // If no occurrences are found, initialize matrices uniformly. // Convert spacer info to probabilities. // Throw out unused motifs.