nodes: 0, module; 1, function_definition; 2, function_name:best_kmers; 3, parameters; 4, block; 5, identifier:dt; 6, identifier:response; 7, identifier:sequence; 8, default_parameter; 9, default_parameter; 10, default_parameter; 11, default_parameter; 12, default_parameter; 13, expression_statement; 14, expression_statement; 15, expression_statement; 16, if_statement; 17, expression_statement; 18, expression_statement; 19, expression_statement; 20, expression_statement; 21, comment:# which coefficients are nonzero?=; 22, expression_statement; 23, comment:# perform stepwise selection; 24, comment:#; 25, comment:# TODO - how do we deal with the intercept?; 26, comment:# largest number of motifs where they don't differ by more than 1 k-mer; 27, function_definition; 28, expression_statement; 29, return_statement; 30, identifier:k; 31, integer:6; 32, identifier:consider_shift; 33, True; 34, identifier:n_cores; 35, integer:1; 36, identifier:seq_align; 37, string:"start"; 38, identifier:trim_seq_len; 39, None; 40, comment:"""
      Find best k-mers for CONCISE initialization.
  
      Args:
          dt (pd.DataFrame): Table containing response variable and sequence.
          response (str): Name of the column used as the reponse variable.
          sequence (str): Name of the column storing the DNA/RNA sequences.
          k (int): Desired k-mer length.
          n_cores (int): Number of cores to use for computation. It can use up to 3 cores.
          consider_shift (boolean): When performing stepwise k-mer selection. Is TATTTA similar to ATTTAG?
          seq_align (str): one of ``{"start", "end"}``. To which end should we align sequences?
          trim_seq_len (int): Consider only first `trim_seq_len` bases of each sequence when generating the sequence design matrix. If :python:`None`, set :py:attr:`trim_seq_len` to the longest sequence length, hence whole sequences are considered.
  
      Returns:
          string list: Best set of motifs for this dataset sorted with respect to
                       confidence (best candidate occuring first).
  
  
      Details:
          First a lasso model gets fitted to get a set of initial motifs. Next, the best
          subset of unrelated motifs is selected by stepwise selection.
      """; 41, assignment; 42, assignment; 43, comparison_operator:trim_seq_len is not None; 44, block; 45, assignment; 46, assignment; 47, assignment; 48, call; 49, assignment; 50, function_name:find_next_best; 51, parameters; 52, block; 53, assignment; 54, identifier:selected_kmers; 55, identifier:y; 56, subscript; 57, identifier:seq; 58, subscript; 59, identifier:trim_seq_len; 60, None; 61, expression_statement; 62, expression_statement; 63, identifier:dt_kmer; 64, call; 65, identifier:Xsp; 66, call; 67, identifier:en; 68, call; 69, attribute; 70, argument_list; 71, identifier:nonzero_kmers; 72, call; 73, identifier:dt_kmer; 74, identifier:y; 75, identifier:selected_kmers; 76, identifier:to_be_selected_kmers; 77, default_parameter; 78, expression_statement; 79, expression_statement; 80, expression_statement; 81, expression_statement; 82, function_definition; 83, expression_statement; 84, if_statement; 85, identifier:selected_kmers; 86, call; 87, identifier:dt; 88, identifier:response; 89, identifier:dt; 90, identifier:sequence; 91, assignment; 92, assignment; 93, identifier:kmer_count; 94, argument_list; 95, identifier:csc_matrix; 96, argument_list; 97, identifier:ElasticNet; 98, argument_list; 99, identifier:en; 100, identifier:fit; 101, identifier:Xsp; 102, identifier:y; 103, attribute; 104, argument_list; 105, identifier:consider_shift; 106, True; 107, comment:"""
          perform stepwise model selection while preventing to add a motif similar to the
          already selected motifs.
          """; 108, assignment; 109, assignment; 110, call; 111, function_name:select_criterion; 112, parameters; 113, block; 114, assignment; 115, comparison_operator:len(to_be_selected_kmers) == 0; 116, block; 117, else_clause; 118, identifier:find_next_best; 119, argument_list; 120, identifier:seq; 121, call; 122, identifier:seq; 123, list_comprehension; 124, identifier:seq; 125, identifier:k; 126, identifier:dt_kmer; 127, keyword_argument; 128, keyword_argument; 129, keyword_argument; 130, subscript; 131, identifier:tolist; 132, pattern_list; 133, call; 134, identifier:kmer; 135, call; 136, attribute; 137, argument_list; 138, identifier:s1; 139, identifier:s2; 140, default_parameter; 141, if_statement; 142, if_statement; 143, if_statement; 144, return_statement; 145, identifier:to_be_selected_kmers; 146, list_comprehension; 147, call; 148, integer:0; 149, return_statement; 150, comment:# regress out the new feature; 151, block; 152, identifier:dt_kmer; 153, identifier:y; 154, list; 155, identifier:nonzero_kmers; 156, identifier:consider_shift; 157, identifier:pad_sequences; 158, argument_list; 159, call; 160, for_in_clause; 161, identifier:alpha; 162, integer:1; 163, identifier:standardize; 164, False; 165, identifier:n_splits; 166, integer:3; 167, attribute; 168, comparison_operator:en.coef_ != 0; 169, identifier:F; 170, identifier:pval; 171, identifier:f_regression; 172, argument_list; 173, attribute; 174, argument_list; 175, identifier:selected_kmers; 176, identifier:append; 177, identifier:kmer; 178, identifier:consider_shift; 179, True; 180, comparison_operator:hamming_distance(s1, s2) <= 1; 181, block; 182, boolean_operator; 183, block; 184, boolean_operator; 185, block; 186, True; 187, identifier:ckmer; 188, for_in_clause; 189, if_clause; 190, identifier:len; 191, argument_list; 192, identifier:selected_kmers; 193, expression_statement; 194, expression_statement; 195, expression_statement; 196, return_statement; 197, identifier:seq; 198, keyword_argument; 199, keyword_argument; 200, attribute; 201, argument_list; 202, identifier:s; 203, identifier:seq; 204, attribute; 205, identifier:values; 206, attribute; 207, integer:0; 208, subscript; 209, identifier:y; 210, identifier:to_be_selected_kmers; 211, identifier:pop; 212, call; 213, call; 214, integer:1; 215, return_statement; 216, identifier:consider_shift; 217, comparison_operator:hamming_distance(s1[1:], s2[:-1]) == 0; 218, return_statement; 219, identifier:consider_shift; 220, comparison_operator:hamming_distance(s1[:-1], s2[1:]) == 0; 221, return_statement; 222, identifier:ckmer; 223, identifier:to_be_selected_kmers; 224, call; 225, identifier:to_be_selected_kmers; 226, assignment; 227, call; 228, assignment; 229, call; 230, identifier:align; 231, identifier:seq_align; 232, identifier:maxlen; 233, identifier:trim_seq_len; 234, identifier:s; 235, identifier:replace; 236, string:"N"; 237, string:""; 238, identifier:dt_kmer; 239, identifier:columns; 240, identifier:en; 241, identifier:coef_; 242, identifier:dt_kmer; 243, identifier:to_be_selected_kmers; 244, attribute; 245, argument_list; 246, identifier:hamming_distance; 247, argument_list; 248, False; 249, call; 250, integer:0; 251, False; 252, call; 253, integer:0; 254, False; 255, identifier:select_criterion; 256, argument_list; 257, identifier:lm; 258, call; 259, attribute; 260, argument_list; 261, identifier:y_new; 262, binary_operator:y - lm.predict(dt_kmer[selected_kmers]); 263, identifier:find_next_best; 264, argument_list; 265, identifier:pval; 266, identifier:argmin; 267, identifier:s1; 268, identifier:s2; 269, identifier:hamming_distance; 270, argument_list; 271, identifier:hamming_distance; 272, argument_list; 273, identifier:ckmer; 274, identifier:kmer; 275, identifier:consider_shift; 276, identifier:LinearRegression; 277, argument_list; 278, identifier:lm; 279, identifier:fit; 280, subscript; 281, identifier:y; 282, identifier:y; 283, call; 284, identifier:dt_kmer; 285, identifier:y_new; 286, identifier:selected_kmers; 287, identifier:to_be_selected_kmers; 288, identifier:consider_shift; 289, subscript; 290, subscript; 291, subscript; 292, subscript; 293, identifier:dt_kmer; 294, identifier:selected_kmers; 295, attribute; 296, argument_list; 297, identifier:s1; 298, slice; 299, identifier:s2; 300, slice; 301, identifier:s1; 302, slice; 303, identifier:s2; 304, slice; 305, identifier:lm; 306, identifier:predict; 307, subscript; 308, integer:1; 309, unary_operator; 310, unary_operator; 311, integer:1; 312, identifier:dt_kmer; 313, identifier:selected_kmers; 314, integer:1; 315, integer:1
code: def best_kmers(dt, response, sequence, k=6, consider_shift=True, n_cores=1,
                 seq_align="start", trim_seq_len=None):
      """
      Find best k-mers for CONCISE initialization.
  
      Args:
          dt (pd.DataFrame): Table containing response variable and sequence.
          response (str): Name of the column used as the reponse variable.
          sequence (str): Name of the column storing the DNA/RNA sequences.
          k (int): Desired k-mer length.
          n_cores (int): Number of cores to use for computation. It can use up to 3 cores.
          consider_shift (boolean): When performing stepwise k-mer selection. Is TATTTA similar to ATTTAG?
          seq_align (str): one of ``{"start", "end"}``. To which end should we align sequences?
          trim_seq_len (int): Consider only first `trim_seq_len` bases of each sequence when generating the sequence design matrix. If :python:`None`, set :py:attr:`trim_seq_len` to the longest sequence length, hence whole sequences are considered.
  
      Returns:
          string list: Best set of motifs for this dataset sorted with respect to
                       confidence (best candidate occuring first).
  
  
      Details:
          First a lasso model gets fitted to get a set of initial motifs. Next, the best
          subset of unrelated motifs is selected by stepwise selection.
      """
      y = dt[response]
      seq = dt[sequence]
  
      if trim_seq_len is not None:
          seq = pad_sequences(seq, align=seq_align, maxlen=trim_seq_len)
          seq = [s.replace("N", "") for s in seq]
  
      dt_kmer = kmer_count(seq, k)
      Xsp = csc_matrix(dt_kmer)
      en = ElasticNet(alpha=1, standardize=False, n_splits=3)
      en.fit(Xsp, y)
      # which coefficients are nonzero?=
      nonzero_kmers = dt_kmer.columns.values[en.coef_ != 0].tolist()
  
      # perform stepwise selection
      #
      # TODO - how do we deal with the intercept?
  
      # largest number of motifs where they don't differ by more than 1 k-mer
      def find_next_best(dt_kmer, y, selected_kmers, to_be_selected_kmers, consider_shift=True):
          """
          perform stepwise model selection while preventing to add a motif similar to the
          already selected motifs.
          """
          F, pval = f_regression(dt_kmer[to_be_selected_kmers], y)
          kmer = to_be_selected_kmers.pop(pval.argmin())
          selected_kmers.append(kmer)
  
          def select_criterion(s1, s2, consider_shift=True):
              if hamming_distance(s1, s2) <= 1:
                  return False
              if consider_shift and hamming_distance(s1[1:], s2[:-1]) == 0:
                  return False
              if consider_shift and hamming_distance(s1[:-1], s2[1:]) == 0:
                  return False
              return True
  
          to_be_selected_kmers = [ckmer for ckmer in to_be_selected_kmers
                                  if select_criterion(ckmer, kmer, consider_shift)]
          if len(to_be_selected_kmers) == 0:
              return selected_kmers
          else:
              # regress out the new feature
              lm = LinearRegression()
              lm.fit(dt_kmer[selected_kmers], y)
              y_new = y - lm.predict(dt_kmer[selected_kmers])
              return find_next_best(dt_kmer, y_new, selected_kmers, to_be_selected_kmers, consider_shift)
  
      selected_kmers = find_next_best(dt_kmer, y, [], nonzero_kmers, consider_shift)
  
      return selected_kmers