Pre Module API Reference

Module for pre-processing to generate LandscapeFiles from ST data.

Network

Bases: object

Clustergrammer.py takes a matrix as input (either from a file of a Pandas DataFrame), normalizes/filters, hierarchically clusters, and produces the :ref:visualization_json for :ref:clustergrammer_js.

Networks have two states:

1. the data state, where they are stored as a matrix and nodes
2. the viz state where they are stored as viz.links, viz.row_nodes, and viz.col_nodes.

The goal is to start in a data-state and produce a viz-state of the network that will be used as input to clustergram.js.

Source code in src/celldega/clust/__init__.py

class Network(object):
  '''
  Clustergrammer.py takes a matrix as input (either from a file of a Pandas DataFrame), normalizes/filters, hierarchically clusters, and produces the :ref:`visualization_json` for :ref:`clustergrammer_js`.

  Networks have two states:

    1. the data state, where they are stored as a matrix and nodes
    2. the viz state where they are stored as viz.links, viz.row_nodes, and viz.col_nodes.

  The goal is to start in a data-state and produce a viz-state of
  the network that will be used as input to clustergram.js.
  '''

  def __init__(self, widget=None):
    initialize_net.main(self, widget)

  def reset(self):
    '''
    This re-initializes the Network object.
    '''
    initialize_net.main(self)

  def load_file(self, filename):
    '''
    Load TSV file.
    '''
    load_data.load_file(self, filename)

  def load_file_as_string(self, file_string, filename=''):
    '''
    Load file as a string.
    '''
    load_data.load_file_as_string(self, file_string, filename=filename)


  def load_stdin(self):
    '''
    Load stdin TSV-formatted string.
    '''
    load_data.load_stdin(self)

  def load_tsv_to_net(self, file_buffer, filename=None):
    '''
    This will load a TSV matrix file buffer; this is exposed so that it will
    be possible to load data without having to read from a file.
    '''
    load_data.load_tsv_to_net(self, file_buffer, filename)

  def load_vect_post_to_net(self, vect_post):
    '''
    Load data in the vector format JSON.
    '''
    load_vect_post.main(self, vect_post)

  def load_data_file_to_net(self, filename):
    '''
    Load Clustergrammer's dat format (saved as JSON).
    '''
    inst_dat = self.load_json_to_dict(filename)
    load_data.load_data_to_net(self, inst_dat)

  def cluster(self, dist_type='cosine', run_clustering=True,
                 dendro=True, views=[],
                 linkage_type='average', sim_mat=False, filter_sim=0.0,
                 calc_cat_pval=False, run_enrichr=None, enrichrgram=None,
                 clust_library='scipy', min_samples=1, min_cluster_size=2):
    '''
    The main function performs hierarchical clustering, optionally generates
    filtered views (e.g. row-filtered views), and generates the :
    ``visualization_json``.

    Used to set views equal to ['N_row_sum', 'N_row_var']
    '''
    initialize_net.viz(self)

    make_clust_fun.make_clust(self, dist_type=dist_type,
                                    run_clustering=run_clustering,
                                    dendro=dendro,
                                    requested_views=views,
                                    linkage_type=linkage_type,
                                    sim_mat=sim_mat,
                                    filter_sim=filter_sim,
                                    calc_cat_pval=calc_cat_pval,
                                    run_enrichr=run_enrichr,
                                    enrichrgram=enrichrgram,
                                    clust_library=clust_library,
                                    min_samples=min_samples,
                                    min_cluster_size=min_cluster_size)

  def swap_nan_for_zero(self):
    '''
    Swaps all NaN (numpy NaN) instances for zero.
    '''
    self.dat['mat'][np.isnan(self.dat['mat'])] = 0

  def load_df(self, df_ini, meta_col=None, meta_row=None, col_cats=None,
              row_cats=None, is_downsampled=False, meta_ds_row=None,
              meta_ds_col=None):
    '''
    Load Pandas DataFrame.
    '''
    self.reset()

    # load dataframe
    df = deepcopy(df_ini)


    if is_downsampled:
      if meta_ds_col is not None:
        self.meta_ds_col = meta_ds_col
      if meta_ds_row is not None:
        self.meta_ds_row = meta_ds_row

    # define downsampled status
    self.is_downsampled = is_downsampled
    # print('load_df: is_downsampled', is_downsampled)

    if hasattr(self, 'meta_col') == False and hasattr(self, 'meta_row') == False:
      self.meta_cat = False

    # load metadata
    if isinstance(meta_col, pd.DataFrame):
      self.meta_col = meta_col

      if col_cats is None:
        self.col_cats = meta_col.columns.tolist()
      else:
        self.col_cats = col_cats

      self.meta_cat = True

    if isinstance(meta_row, pd.DataFrame):
      self.meta_row = meta_row

      if row_cats is None:
        self.row_cats = meta_row.columns.tolist()
      else:
        self.row_cats = row_cats

      self.meta_cat = True

    data_formats.df_to_dat(self, df, define_cat_colors=True)

  def export_df(self):
    '''
    Export Pandas DataFrame/
    '''
    df = data_formats.dat_to_df(self)

    # drop tuple categories if downsampling
    if self.is_downsampled:
      df.columns = self.dat['nodes']['col']
      df.index = self.dat['nodes']['row']

    return df

  def df_to_dat(self, df, define_cat_colors=False):
    '''
    Load Pandas DataFrame (will be deprecated).
    '''
    data_formats.df_to_dat(self, df, define_cat_colors)

  def set_matrix_colors(self, pos='red', neg='blue'):

    self.viz['matrix_colors'] = {}
    self.viz['matrix_colors']['pos'] = pos
    self.viz['matrix_colors']['neg'] = neg

  def set_global_cat_colors(self, df_meta):

    for inst_name in df_meta.index.tolist():
      inst_color = df_meta.loc[inst_name, 'color']

      self.viz['global_cat_colors'][inst_name] = inst_color

  def set_cat_color(self, axis, cat_index, cat_name, inst_color):

    if axis == 0:
      axis = 'row'
    if axis == 1:
      axis = 'col'

    try:
      # process cat_index
      cat_index = cat_index - 1
      cat_index = 'cat-' + str(cat_index)

      self.viz['cat_colors'][axis][cat_index][cat_name] = inst_color

    except:
      print('there was an error setting the category color')

  def dat_to_df(self):
    '''
    Export Pandas DataFrams (will be deprecated).
    '''
    return data_formats.dat_to_df(self)

  def export_net_json(self, net_type='viz', indent='no-indent'):
    '''
    Export dat or viz JSON.
    '''
    return export_data.export_net_json(self, net_type, indent)

  def export_viz_to_widget(self, which_viz='viz'):
    '''
    Export viz JSON, for use with clustergrammer_widget. Formerly method was
    named widget.
    '''

    return export_data.export_net_json(self, which_viz, 'no-indent')

  def widget(self, which_viz='viz', link_net=None, link_net_js=None, clust_library='scipy',
    min_samples=1, min_cluster_size=2):
    '''
    Generate a widget visualization using the widget. The export_viz_to_widget
    method passes the visualization JSON to the instantiated widget, which is
    returned and visualized on the front-end.
    '''
    # run clustering if necessary
    if len(self.viz['row_nodes']) == 0:
      self.cluster(clust_library=clust_library, min_samples=min_samples,
                   min_cluster_size=min_cluster_size)

      # add manual_category to viz json
      if 'manual_category' in self.dat:
        self.viz['manual_category'] = self.dat['manual_category']

      # add pre-z-score data to viz
      if 'pre_zscore' in self.dat:
        self.viz['pre_zscore'] = self.dat['pre_zscore']

    self.widget_instance = self.widget_class(network = self.export_viz_to_widget(which_viz))

    # initialize manual category
    if 'manual_category' in self.dat:
      manual_cat = {}
      axis = 'col'
      manual_cat[axis] = {}
      manual_cat[axis]['col_cat_colors'] = self.viz['cat_colors'][axis]['cat-0']

      man_cat_name = self.dat['manual_category'][axis]
      if axis == 'col':
        manual_cat[axis][man_cat_name] = self.meta_col[man_cat_name].to_dict()

      self.widget_instance.manual_cat = json.dumps(manual_cat)

      self.widget_instance.observe(self.get_manual_category, names='manual_cat')

    # add link (python)
    if link_net is not None:
      inst_link = widgets.link(
                                (self.widget_instance, 'manual_cat'),
                                (link_net.widget_instance, 'manual_cat')
                               )
      self.widget_instance.link = inst_link

    # add jslink (JavaScript)
    if link_net_js is not None:

      inst_link = widgets.jslink(
                                (self.widget_instance, 'manual_cat'),
                                (link_net_js.widget_instance, 'manual_cat')
                               )
      self.widget_instance.link = inst_link

    return self.widget_instance


  def widget_df(self):
    '''
    Export a DataFrame from the front-end visualization. For instance, a user
    can filter to show only a single cluster using the dendrogram and then
    get a dataframe of this cluster using the widget_df method.
    '''

    if hasattr(self, 'widget_instance') == True:

      if self.widget_instance.mat_string != '':

        tmp_net = deepcopy(Network())

        df_string = self.widget_instance.mat_string

        tmp_net.load_file_as_string(df_string)

        df = tmp_net.export_df()

        return df

      else:
        return self.export_df()

    else:
      if hasattr(self, 'widget_class') == True:
        print('Please make the widget before exporting the widget DataFrame.')
        print('Do this using the widget method: net.widget()')

      else:
        print('Can not make widget because Network has no attribute widget_class')
        print('Please instantiate Network with clustergrammer_widget using: Network(clustergrammer_widget)')

  def write_json_to_file(self, net_type, filename, indent='no-indent'):
    '''
    Save dat or viz as a JSON to file.
    '''
    export_data.write_json_to_file(self, net_type, filename, indent)

  def write_matrix_to_tsv(self, filename=None, df=None):
    '''
    Export data-matrix to file.
    '''
    return export_data.write_matrix_to_tsv(self, filename, df)

  def filter_sum(self, threshold, take_abs=True, axis=None, inst_rc=None):
    '''
    Filter a network's rows or columns based on the sum across rows or columns.
    '''

    if axis is None:
      axis = inst_rc
      print('warning inst_rc argument will be deprecated, please use axis')

      if inst_rc is None:
        print('please provide axis argument')

    inst_df = self.dat_to_df()
    if axis == 'row':
      inst_df = run_filter.df_filter_row_sum(inst_df, threshold, take_abs)
    elif axis == 'col':
      inst_df = run_filter.df_filter_col_sum(inst_df, threshold, take_abs)
    self.df_to_dat(inst_df)

  def filter_N_top(self, N_top, rank_type='sum', inst_rc=None, axis=None):
    '''
    Filter the matrix rows or columns based on sum/variance, and only keep the top
    N.
    '''

    if axis is None:
      axis = inst_rc
      print('warning inst_rc argument will be deprecated, please use axis')

      if inst_rc is None:
        print('please provide axis argument')

    inst_df = self.dat_to_df()

    inst_df = run_filter.filter_N_top(inst_rc, inst_df, N_top, rank_type)

    self.df_to_dat(inst_df)

  def filter_threshold(self, threshold, num_occur=1, inst_rc=None, axis=None):
    '''
    Filter the matrix rows or columns based on num_occur values being above a
    threshold (in absolute value).
    '''

    if axis is None:
      axis = inst_rc
      print('warning inst_rc argument will be deprecated, please use axis')

      if inst_rc is None:
        print('please provide axis argument')

    inst_df = self.dat_to_df()

    inst_df = run_filter.filter_threshold(inst_df, axis, threshold,
      num_occur)

    self.df_to_dat(inst_df)

  def filter_cat(self, axis, cat_index, cat_name):
    '''
    Filter the matrix based on their category. cat_index is the index of the category, the first category has index=1.
    '''
    run_filter.filter_cat(self, axis, cat_index, cat_name)

  def filter_names(self, axis, names):
    '''
    Filter the visualization using row/column names. The function takes, axis ('row'/'col') and names, a list of strings.
    '''
    run_filter.filter_names(self, axis, names)

  def clip(self, lower=None, upper=None):
    '''
    Trim values at input thresholds using pandas function
    '''
    df = self.export_df()
    df = df.clip(lower=lower, upper=upper)
    self.load_df(df)

  def normalize(self, df=None, norm_type='zscore', axis='row', z_clip=None):
    '''
    Normalize the matrix rows or columns using Z-score (zscore) or Quantile Normalization (qn). Users can optionally pass in a DataFrame to be normalized (and this will be incorporated into the Network object).
    '''
    normalize_fun.run_norm(self, df, norm_type, axis, z_clip)

  def downsample(self, df=None, ds_type='kmeans', axis='row', num_samples=100,
                 random_state=1000, ds_name='Downsample',
                 ds_cluster_name='cluster'):
    '''
    Downsample the matrix rows or columns (currently supporting kmeans only). Users can optionally pass in a DataFrame to be downsampled (and this will be incorporated into the network object).
    '''
    return downsample_fun.main(self, df, ds_type, axis, num_samples, random_state, ds_name, ds_cluster_name)

  def random_sample(self, num_samples, df=None, replace=False, weights=None, random_state=100, axis='row'):
    '''
    Return random sample of matrix.
    '''

    if df is None:
      df = self.dat_to_df()

    if axis == 'row':
      axis = 0
    if axis == 'col':
      axis = 1

    df = self.export_df()
    df = df.sample(n=num_samples, replace=replace, weights=weights, random_state=random_state,  axis=axis)

    self.load_df(df)

  def add_cats(self, axis, cat_data):
    '''
    Add categories to rows or columns using cat_data array of objects. Each object in cat_data is a dictionary with one key (category title) and value (rows/column names) that have this category. Categories will be added onto the existing categories and will be added in the order of the objects in the array.

    Example ``cat_data``::


        [
          {
            "title": "First Category",
            "cats": {
              "true": [
                "ROS1",
                "AAK1"
              ]
            }
          },
          {
            "title": "Second Category",
            "cats": {
              "something": [
                "PDK4"
              ]
            }
          }
        ]


    '''
    for inst_data in cat_data:
      categories.add_cats(self, axis, inst_data)

  def Iframe_web_app(self, filename=None, width=1000, height=800):

    link = iframe_web_app.main(self, filename, width, height)

    return link

  def enrichrgram(self, lib, axis='row'):
    '''
    Add Enrichr gene enrichment results to your visualization (where your rows
    are genes). Run enrichrgram before clustering to incldue enrichment results
    as row categories. Enrichrgram can also be run on the front-end using the
    Enrichr logo at the top left.

    Set lib to the Enrichr library that you want to use for enrichment analysis.
    Libraries included:

      * ChEA_2016
      * KEA_2015
      * ENCODE_TF_ChIP-seq_2015
      * ENCODE_Histone_Modifications_2015
      * Disease_Perturbations_from_GEO_up
      * Disease_Perturbations_from_GEO_down
      * GO_Molecular_Function_2015
      * GO_Biological_Process_2015
      * GO_Cellular_Component_2015
      * Reactome_2016
      * KEGG_2016
      * MGI_Mammalian_Phenotype_Level_4
      * LINCS_L1000_Chem_Pert_up
      * LINCS_L1000_Chem_Pert_down

    '''

    df = self.export_df()
    df, bar_info = enr_fun.add_enrichr_cats(df, axis, lib)
    self.load_df(df)

    self.dat['enrichrgram_lib'] = lib
    self.dat['row_cat_bars'] = bar_info

  @staticmethod
  def load_gmt(filename):
    return load_data.load_gmt(filename)

  @staticmethod
  def load_json_to_dict(filename):
    return load_data.load_json_to_dict(filename)

  @staticmethod
  def save_dict_to_json(inst_dict, filename, indent='no-indent'):
    export_data.save_dict_to_json(inst_dict, filename, indent)

  @staticmethod
  def load_gene_exp_to_df(inst_path):
    '''
    Loads gene expression data from 10x in sparse matrix format and returns a
    Pandas dataframe
    '''

    import pandas as pd
    from scipy import io
    from scipy import sparse
    from ast import literal_eval as make_tuple

    # matrix
    Matrix = io.mmread( inst_path + 'matrix.mtx')
    mat = Matrix.todense()

    # genes
    filename = inst_path + 'genes.tsv'
    f = open(filename, 'r')
    lines = f.readlines()
    f.close()

    # # add unique id to all genes
    # genes = []
    # unique_id = 0
    # for inst_line in lines:
    #     inst_line = inst_line.strip().split()

    #     if len(inst_line) > 1:
    #       inst_gene = inst_line[1]
    #     else:
    #       inst_gene = inst_line[0]

    #     genes.append(inst_gene + '_' + str(unique_id))
    #     unique_id = unique_id + 1

    # add unique id only to duplicate genes
    ini_genes = []
    for inst_line in lines:
        inst_line = inst_line.strip().split()
        if len(inst_line) > 1:
          inst_gene = inst_line[1]
        else:
          inst_gene = inst_line[0]
        ini_genes.append(inst_gene)

    gene_name_count = pd.Series(ini_genes).value_counts()
    duplicate_genes = gene_name_count[gene_name_count > 1].index.tolist()

    dup_index = {}
    genes = []
    for inst_row in ini_genes:

      # add index to non-unique genes
      if inst_row in duplicate_genes:

        # calc_non-unque index
        if inst_row not in dup_index:
          dup_index[inst_row] = 1
        else:
          dup_index[inst_row] = dup_index[inst_row] + 1

        new_row = inst_row + '_' + str(dup_index[inst_row])

      else:
        new_row = inst_row

      genes.append(new_row)

    # barcodes
    filename = inst_path + 'barcodes.tsv'
    f = open(filename, 'r')
    lines = f.readlines()
    f.close()

    cell_barcodes = []
    for inst_bc in lines:
        inst_bc = inst_bc.strip().split('\t')

        # remove dash from barcodes if necessary
        if '-' in inst_bc[0]:
          inst_bc[0] = inst_bc[0].split('-')[0]

        cell_barcodes.append(inst_bc[0])

    # parse tuples if necessary
    try:
        cell_barcodes = [make_tuple(x) for x in cell_barcodes]
    except:
        pass

    try:
        genes = [make_tuple(x) for x in genes]
    except:
        pass

    # make dataframe
    df = pd.DataFrame(mat, index=genes, columns=cell_barcodes)

    return df

  @staticmethod
  def save_gene_exp_to_mtx_dir(inst_path, df):

    import os
    from scipy import io
    from scipy import sparse

    if not os.path.exists(inst_path):
        os.makedirs(inst_path)

    genes = df.index.tolist()
    barcodes = df.columns.tolist()

    save_list_to_tsv(genes, inst_path + 'genes.tsv')
    save_list_to_tsv(barcodes, inst_path + 'barcodes.tsv')

    mat_ge = df.values
    mat_ge_sparse = sparse.coo_matrix(mat_ge)

    io.mmwrite( inst_path + 'matrix.mtx', mat_ge_sparse)

  @staticmethod
  def save_list_to_tsv(inst_list, filename):

      f = open(filename, 'w')
      for inst_line in inst_list:
          f.write(str(inst_line) + '\n')

      f.close()

  def sim_same_and_diff_category_samples(self, df, cat_index=1, dist_type='cosine',
                                         equal_var=False, plot_roc=True,
                                         precalc_dist=False, calc_roc=True):
      '''
      Calculate the similarity of samples from the same and different categories. The
      cat_index gives the index of the category, where 1 in the first category
      '''

      cols = df.columns.tolist()

      if type(precalc_dist) == bool:
          # compute distnace between rows (transpose to get cols as rows)
          dist_arr = 1 - pdist(df.transpose(), metric=dist_type)
      else:
          dist_arr = precalc_dist

      # generate sample names with categories
      sample_combos = list(combinations(range(df.shape[1]),2))

      sample_names = [str(ind) + '_same' if cols[x[0]][cat_index] == cols[x[1]][cat_index] else str(ind) + '_different' for ind, x in enumerate(sample_combos)]

      ser_dist = pd.Series(data=dist_arr, index=sample_names)

      # find same-cat sample comparisons
      same_cat = [x for x in sample_names if x.split('_')[1] == 'same']

      # find diff-cat sample comparisons
      diff_cat = [x for x in sample_names if x.split('_')[1] == 'different']

      # make series of same and diff category sample comparisons
      ser_same = ser_dist[same_cat]
      ser_same.name = 'Same Category'
      ser_diff = ser_dist[diff_cat]
      ser_diff.name = 'Different Category'

      sim_dict = {}
      roc_data = {}
      sim_data = {}

      sim_dict['same'] = ser_same
      sim_dict['diff'] = ser_diff

      pval_dict = {}
      ttest_stat, pval_dict['ttest'] = ttest_ind(ser_diff, ser_same, equal_var=equal_var)

      ttest_stat, pval_dict['mannwhitney'] = mannwhitneyu(ser_diff, ser_same)

      if calc_roc:
          # calc AUC
          true_index = list(np.ones(sim_dict['same'].shape[0]))
          false_index = list(np.zeros(sim_dict['diff'].shape[0]))
          y_true = true_index + false_index

          true_val = list(sim_dict['same'].values)
          false_val = list(sim_dict['diff'].values)
          y_score = true_val + false_val

          fpr, tpr, thresholds = roc_curve(y_true, y_score)

          inst_auc = auc(fpr, tpr)

          if plot_roc:
              plt.figure()
              plt.plot(fpr, tpr)
              plt.plot([0, 1], [0, 1], color='navy', linestyle='--')
              plt.figure(figsize=(10,10))

              print('AUC', inst_auc)

          roc_data['true'] = y_true
          roc_data['score'] = y_score
          roc_data['fpr'] = fpr
          roc_data['tpr'] = tpr
          roc_data['thresholds'] = thresholds
          roc_data['auc'] = inst_auc

      sim_data['sim_dict'] = sim_dict
      sim_data['pval_dict'] = pval_dict
      sim_data['roc_data'] = roc_data

      return sim_data

  def generate_signatures(self, df_data, df_meta, category_name, pval_cutoff=0.05,
                          num_top_dims=False, verbose=True, equal_var=False):

      '''
      Generate signatures for column categories
      T-test is run for each category in a one-vs-all manner. The num_top_dims overrides
      the P-value cutoff.
      '''

      df_t = df_data.transpose()

      # remove columns (dimensions) with constant values
      orig_num_cols = df_t.shape[1]
      df_t = df_t.loc[:, (df_t != df_t.iloc[0]).any()]
      if df_t.shape[1] < orig_num_cols:
        print('dropped columns with constant values')

      # make tuple rows
      df_t.index = [(x, df_meta.loc[x, category_name]) for x in df_t.index.tolist()]
      category_level = 1

      df = self.row_tuple_to_multiindex(df_t)

      cell_types = sorted(list(set(df.index.get_level_values(category_level).tolist())))

      keep_genes = []
      keep_genes_dict = {}
      gene_pval_dict = {}
      all_fold_info = {}

      for inst_ct in cell_types:

          inst_ct_mat = df.xs(key=inst_ct, level=category_level)
          inst_other_mat = df.drop(inst_ct, level=category_level)

          # save mean values and fold change
          fold_info = {}
          fold_info['cluster_mean'] = inst_ct_mat.mean()
          fold_info['other_mean'] = inst_other_mat.mean()
          fold_info['log2_fold'] = fold_info['cluster_mean']/fold_info['other_mean']
          fold_info['log2_fold'] = fold_info['log2_fold'].apply(np.log2)
          all_fold_info[inst_ct] = fold_info

          inst_stats, inst_pvals = ttest_ind(inst_ct_mat, inst_other_mat, axis=0, equal_var=equal_var)

          ser_pval = pd.Series(data=inst_pvals, index=df.columns.tolist()).sort_values()

          if num_top_dims == False:
              ser_pval_keep = ser_pval[ser_pval < pval_cutoff]
          else:
              ser_pval_keep = ser_pval[:num_top_dims]

          gene_pval_dict[inst_ct] = ser_pval_keep

          inst_keep = ser_pval_keep.index.tolist()
          keep_genes.extend(inst_keep)
          keep_genes_dict[inst_ct] = inst_keep

      keep_genes = sorted(list(set(keep_genes)))

      df_gbm = df.groupby(level=category_level).mean().transpose()
      cols = df_gbm.columns.tolist()

      df_sig = df_gbm.loc[keep_genes]

      if len(keep_genes) == 0 and verbose:
          print('found no informative dimensions')

      df_gene_pval = pd.concat(gene_pval_dict, axis=1, sort=False)

      df_diff = {}
      for inst_col in df_gene_pval:

          inst_pvals = df_gene_pval[inst_col]

          # nans represent dimensions that did not meet pval or top threshold
          inst_pvals = inst_pvals.dropna().sort_values()

          inst_genes = inst_pvals.index.tolist()

          # prevent failure if no cells have this category
          if inst_pvals.shape[0] > 0:
              rej, pval_corr = smm.multipletests(inst_pvals, 0.05, method='fdr_bh')[:2]

              ser_pval = pd.Series(data=inst_pvals, index=inst_genes, name='P-values')
              ser_bh_pval = pd.Series(data=pval_corr, index=inst_genes, name='BH P-values')
              ser_log2_fold = all_fold_info[inst_col]['log2_fold'].loc[inst_genes]
              ser_log2_fold.name = 'Log2 Fold Change'
              ser_cluster_mean = all_fold_info[inst_col]['cluster_mean'].loc[inst_genes]
              ser_cluster_mean.name = 'Cluster Mean'
              ser_other_mean = all_fold_info[inst_col]['other_mean'].loc[inst_genes]
              ser_other_mean.name = 'All Other Mean'
              inst_df = pd.concat([ser_pval, ser_bh_pval, ser_log2_fold, ser_cluster_mean, ser_other_mean], axis=1)

              df_diff[inst_col] = inst_df

      return df_sig, df_diff

  def old_generate_signatures(self, df_ini, category_level, pval_cutoff=0.05,
                          num_top_dims=False, verbose=True, equal_var=False):

      ''' Generate signatures for column categories '''

      # change this to use category name and set caetgory level to 1 always
      ###################################

      df_t = df_ini.transpose()

      # remove columns with constant values
      df_t = df_t.loc[:, (df_t != df_t.iloc[0]).any()]

      df = self.row_tuple_to_multiindex(df_t)

      cell_types = sorted(list(set(df.index.get_level_values(category_level).tolist())))

      keep_genes = []
      keep_genes_dict = {}
      gene_pval_dict = {}
      all_fold_info = {}

      for inst_ct in cell_types:

          inst_ct_mat = df.xs(key=inst_ct, level=category_level)
          inst_other_mat = df.drop(inst_ct, level=category_level)

          # save mean values and fold change
          fold_info = {}
          fold_info['cluster_mean'] = inst_ct_mat.mean()
          fold_info['other_mean'] = inst_other_mat.mean()
          fold_info['log2_fold'] = fold_info['cluster_mean']/fold_info['other_mean']
          fold_info['log2_fold'] = fold_info['log2_fold'].apply(np.log2)
          all_fold_info[inst_ct] = fold_info

          inst_stats, inst_pvals = ttest_ind(inst_ct_mat, inst_other_mat, axis=0, equal_var=equal_var)

          ser_pval = pd.Series(data=inst_pvals, index=df.columns.tolist()).sort_values()

          if num_top_dims == False:
              ser_pval_keep = ser_pval[ser_pval < pval_cutoff]
          else:
              ser_pval_keep = ser_pval[:num_top_dims]

          gene_pval_dict[inst_ct] = ser_pval_keep

          inst_keep = ser_pval_keep.index.tolist()
          keep_genes.extend(inst_keep)
          keep_genes_dict[inst_ct] = inst_keep

      keep_genes = sorted(list(set(keep_genes)))

      df_gbm = df.groupby(level=category_level).mean().transpose()
      cols = df_gbm.columns.tolist()
      new_cols = []
      for inst_col in cols:
          new_col = (inst_col, category_level + ': ' + inst_col)
          new_cols.append(new_col)
      df_gbm.columns = new_cols

      df_sig = df_gbm.loc[keep_genes]

      if len(keep_genes) == 0 and verbose:
          print('found no informative dimensions')

      df_gene_pval = pd.concat(gene_pval_dict, axis=1, sort=False)

      return df_sig, df_gene_pval, all_fold_info

  def predict_cats_from_sigs(self, df_data_ini, df_meta, df_sig_ini,
                             predict='Predicted Category', dist_type='cosine',
                             unknown_thresh=-1):
      ''' Predict category using signature '''

      keep_rows = df_sig_ini.index.tolist()
      data_rows = df_data_ini.index.tolist()

      common_rows = list(set(data_rows).intersection(keep_rows))

      df_data = deepcopy(df_data_ini.loc[common_rows])
      df_sig = deepcopy(df_sig_ini.loc[common_rows])

      # calculate sim_mat of df_data and df_sig
      cell_types = df_sig.columns.tolist()
      barcodes = df_data.columns.tolist()
      sim_mat = 1 - pairwise_distances(df_sig.transpose(), df_data.transpose(), metric=dist_type)
      df_sim = pd.DataFrame(data=sim_mat, index=cell_types, columns=barcodes).transpose()

      # get the top column value (most similar signature)
      df_sim_top = df_sim.idxmax(axis=1)

      # get the maximum similarity of a cell to a cell type definition
      max_sim = df_sim.max(axis=1)

      unknown_cells = max_sim[max_sim < unknown_thresh].index.tolist()

      # assign unknown cells (need category of same name)
      df_sim_top[unknown_cells] = 'Unknown'

      # add predicted category name to top list
      # top_list = df_sim_top.values
      df_sim = df_sim.transpose()

      df_meta[predict] = df_sim_top

      return df_sim, df_meta


  def old_predict_cats_from_sigs(self, df_data_ini, df_sig_ini, dist_type='cosine', predict_level='Predict Category',
                             truth_level=1, unknown_thresh=-1):
      ''' Predict category using signature '''

      keep_rows = df_sig_ini.index.tolist()
      data_rows = df_data_ini.index.tolist()

      common_rows = list(set(data_rows).intersection(keep_rows))

      df_data = deepcopy(df_data_ini.loc[common_rows])
      df_sig = deepcopy(df_sig_ini.loc[common_rows])

      # calculate sim_mat of df_data and df_sig
      cell_types = df_sig.columns.tolist()
      barcodes = df_data.columns.tolist()
      sim_mat = 1 - pairwise_distances(df_sig.transpose(), df_data.transpose(), metric=dist_type)
      df_sim = pd.DataFrame(data=sim_mat, index=cell_types, columns=barcodes).transpose()

      # get the top column value (most similar signature)
      df_sim_top = df_sim.idxmax(axis=1)

      # get the maximum similarity of a cell to a cell type definition
      max_sim = df_sim.max(axis=1)

      unknown_cells = max_sim[max_sim < unknown_thresh].index.tolist()

      # assign unknown cells (need category of same name)
      df_sim_top[unknown_cells] = 'Unknown'

      # add predicted category name to top list
      top_list = df_sim_top.values
      top_list = [ predict_level + ': ' + x[0] if type(x) is tuple else predict_level + ': ' + x  for x in top_list]

      # add cell type category to input data
      df_cat = deepcopy(df_data)
      cols = df_cat.columns.tolist()
      new_cols = []

      # check whether the columns have the true category available
      has_truth = False
      if type(cols[0]) is tuple:
          has_truth = True

      if has_truth:
          new_cols = [tuple(list(a) + [b]) for a,b in zip(cols, top_list)]
      else:
          new_cols = [tuple([a] + [b]) for a,b in zip(cols, top_list)]

      # transfer new categories
      df_cat.columns = new_cols

      # keep track of true and predicted labels
      y_info = {}
      y_info['true'] = []
      y_info['pred'] = []

      if has_truth:
          y_info['true'] = [x[truth_level].split(': ')[1] for x in cols]
          y_info['pred'] = [x.split(': ')[1] for x in top_list]

      return df_cat, df_sim.transpose(), y_info

  def assess_prediction(self, df_meta, truth, pred):
      ''' Generate confusion matrix from y_info '''

      y_info = {}
      y_info['true'] = df_meta[truth].values.tolist()
      y_info['pred'] = df_meta[pred].values.tolist()

      sorted_cats = sorted(list(set(y_info['true'] + y_info['pred'])))
      conf_mat = confusion_matrix(y_info['true'], y_info['pred'], labels=sorted_cats)

      # true columns and pred rows
      df_conf = pd.DataFrame(conf_mat, index=sorted_cats, columns=sorted_cats).transpose()

      total_correct = np.trace(df_conf)
      total_pred = df_conf.sum().sum()
      fraction_correct = total_correct/float(total_pred)

      # calculate ser_correct
      correct_list = []
      cat_counts = df_conf.sum(axis=0)
      all_cols = df_conf.columns.tolist()
      for inst_cat in all_cols:
          inst_correct = df_conf.loc[inst_cat, inst_cat] / cat_counts[inst_cat]
          correct_list.append(inst_correct)

      ser_correct = pd.Series(data=correct_list, index=all_cols)

      return df_conf, ser_correct, fraction_correct

  def old_confusion_matrix_and_correct_series(self, y_info):
      ''' Generate confusion matrix from y_info '''


      a = deepcopy(y_info['true'])
      true_count = dict((i, a.count(i)) for i in set(a))

      a = deepcopy(y_info['pred'])
      pred_count = dict((i, a.count(i)) for i in set(a))

      sorted_cats = sorted(list(set(y_info['true'] + y_info['pred'])))
      conf_mat = confusion_matrix(y_info['true'], y_info['pred'], sorted_cats)
      df_conf = pd.DataFrame(conf_mat, index=sorted_cats, columns=sorted_cats)

      total_correct = np.trace(df_conf)
      total_pred = df_conf.sum().sum()
      fraction_correct = total_correct/float(total_pred)

      # calculate ser_correct
      correct_list = []
      cat_counts = df_conf.sum(axis=1)
      all_cols = df_conf.columns.tolist()
      for inst_cat in all_cols:
          inst_correct = df_conf[inst_cat].loc[inst_cat] / cat_counts[inst_cat]
          correct_list.append(inst_correct)

      ser_correct = pd.Series(data=correct_list, index=all_cols)

      populations = {}
      populations['true'] = true_count
      populations['pred'] = pred_count

      return df_conf, populations, ser_correct, fraction_correct

  def compare_performance_to_shuffled_labels(self, df_data, category_level, num_shuffles=100,
                                             random_seed=99, pval_cutoff=0.05, dist_type='cosine',
                                             num_top_dims=False, predict_level='Predict Category',
                                             truth_level=1, unknown_thresh=-1, equal_var=False,
                                             performance_type='prediction'):
      random.seed(random_seed)

      perform_list = []
      num_shuffles = num_shuffles

      # pre-calculate the distance matrix (similarity matrix) if necessary
      if performance_type == 'cat_sim_auc':
          dist_arr = 1 - pdist(df_data.transpose(), metric=dist_type)

      for inst_run in range(num_shuffles + 1):

          cols = df_data.columns.tolist()
          rows = df_data.index.tolist()
          mat = df_data.values

          shuffled_cols = deepcopy(cols)
          random.shuffle(shuffled_cols)

          # do not perform shuffling the first time to confirm that we get the same
          # results as the unshuffled dataaset
          if inst_run == 0:
              df_shuffle = deepcopy(df_data)
          else:
              df_shuffle = pd.DataFrame(data=mat, columns=shuffled_cols, index=rows)

          # generate signature on shuffled data
          df_sig, keep_genes, keep_genes_dict, fold_info = generate_signatures(df_shuffle,
                                                             category_level,
                                                             pval_cutoff=pval_cutoff,
                                                             num_top_dims=num_top_dims,
                                                             equal_var=equal_var)

          # predictive performance
          if performance_type == 'prediction':

              # predict categories from signature
              df_pred_cat, df_sig_sim, y_info = self.predict_cats_from_sigs(df_shuffle, df_sig,
                  dist_type=dist_type, predict_level=predict_level, truth_level=truth_level,
                  unknown_thresh=unknown_thresh)

              # calc confusion matrix and performance
              df_conf, populations, ser_correct, fraction_correct = self.confusion_matrix_and_correct_series(y_info)

              # store performances of shuffles
              if inst_run > 0:
                  perform_list.append(fraction_correct)
              else:
                  real_performance = fraction_correct
                  print('performance (fraction correct) of unshuffled: ' + str(fraction_correct))

          elif performance_type == 'cat_sim_auc':

              # predict categories from signature
              sim_dict, pval_dict, roc_data = self.sim_same_and_diff_category_samples(df_shuffle,
                  cat_index=1, plot_roc=False, equal_var=equal_var, precalc_dist=dist_arr)

              # store performances of shuffles
              if inst_run > 0:
                  perform_list.append(roc_data['auc'])
              else:
                  real_performance = roc_data['auc']
                  print('performance (category similarity auc) of unshuffled: ' + str(roc_data['auc']))

      perform_ser = pd.Series(perform_list)

      in_top_fraction = perform_ser[perform_ser > real_performance].shape[0]/num_shuffles
      print('real data performs in the top ' + str(in_top_fraction*100) + '% of shuffled labels\n')

      return perform_ser

  @staticmethod
  def box_scatter_plot(df, group, columns=False, rand_seed=100, alpha=0.5,
      dot_color='red', num_row=None, num_col=1, figsize=(10,10),
      start_title='Variable Measurements Across', end_title='Groups',
      group_list=False):

      from scipy import stats
      import pandas as pd

      import matplotlib.pyplot as plt
      # %matplotlib inline

      if columns == False:
          columns = df.columns.tolist()

      plt.figure(figsize=figsize)
      figure_title = start_title + ' ' + group + ' ' + end_title
      plt.suptitle(figure_title, fontsize=20)

      # list of arranged dataframes
      dfs = {}

      for col_num in range(len(columns)):
          column = columns[col_num]
          plot_id = col_num + 1

          # group by column name or multiIndex name
          if group in df.columns.tolist():
              grouped = df.groupby(group)
          else:
              grouped = df.groupby(level=group)

          names, vals, xs = [], [] ,[]

          if type(column) is tuple:
              column_title = column[0]
          else:
              column_title = column

          for i, (name, subdf) in enumerate(grouped):

              names.append(name)

              inst_ser = subdf[column]

              column_name = column_title + '-' + str(name)

              inst_ser.name = column_name
              vals.append(inst_ser)

              np.random.seed(rand_seed)
              xs.append(np.random.normal(i+1, 0.04, subdf.shape[0]))

          ax = plt.subplot(num_row, num_col, plot_id)

          plt.boxplot(vals, labels=names)

          ngroup = len(vals)

          clevels = np.linspace(0., 1., ngroup)

          for x, val, clevel in zip(xs, vals, clevels):

              plt.subplot(num_row, num_col, plot_id)
              plt.scatter(x, val, c=dot_color, alpha=alpha)


          df_arranged = pd.concat(vals, axis=1)

          # anova
          anova_data = [df_arranged[col].dropna() for col in df_arranged]
          f_val, pval = stats.f_oneway(*anova_data)

          if pval < 0.01:
              ax.set_title(column_title + ' P-val: ' + '{:.2e}'.format(pval))
          else:
              pval = round(pval * 100000)/100000
              ax.set_title(column_title + ' P-val: ' + str(pval))

          dfs[column] = df_arranged

      return dfs

  @staticmethod
  def rank_cols_by_anova_pval(df, group, columns=False, rand_seed=100, alpha=0.5, dot_color='red', num_row=None, num_col=1,
                       figsize=(10,10)):

      from scipy import stats
      import numpy as np
      import pandas as pd

      # import matplotlib.pyplot as plt
      # %matplotlib inline

      if columns == False:
          columns = df.columns.tolist()

      # plt.figure(figsize=figsize)

      # list of arranged dataframes
      dfs = {}

      pval_list = []

      for col_num in range(len(columns)):
          column = columns[col_num]
          plot_id = col_num + 1

          # group by column name or multiIndex name
          if group in df.columns.tolist():
              grouped = df.groupby(group)
          else:
              grouped = df.groupby(level=group)

          names, vals, xs = [], [] ,[]

          if type(column) is tuple:
              column_title = column[0]
          else:
              column_title = column

          for i, (name, subdf) in enumerate(grouped):
              names.append(name)

              inst_ser = subdf[column]

              column_name = column_title + '-' + str(name)

              inst_ser.name = column_name
              vals.append(inst_ser)

              np.random.seed(rand_seed)
              xs.append(np.random.normal(i+1, 0.04, subdf.shape[0]))


          ngroup = len(vals)

          df_arranged = pd.concat(vals, axis=1)

          # anova
          anova_data = [df_arranged[col].dropna() for col in df_arranged]
          f_val, pval = stats.f_oneway(*anova_data)

          pval_list.append(pval)

      pval_ser = pd.Series(data=pval_list, index=columns)
      pval_ser = pval_ser.sort_values(ascending=True)

      return pval_ser


  def row_tuple_to_multiindex(self, df):

      import pandas as pd

      from copy import deepcopy
      df_mi = deepcopy(df)
      rows = df_mi.index.tolist()
      titles = []
      for inst_part in rows[0]:

          if ': ' in inst_part:
              inst_title = inst_part.split(': ')[0]
          else:
              inst_title = 'Name'
          titles.append(inst_title)

      new_rows = []
      for inst_row in rows:
          inst_row = list(inst_row)
          new_row = []
          for inst_part in inst_row:
              if ': ' in inst_part:
                  inst_part = inst_part.split(': ')[1]
              new_row.append(inst_part)
          new_row = tuple(new_row)
          new_rows.append(new_row)

      df_mi.index = new_rows

      df_mi.index = pd.MultiIndex.from_tuples(df_mi.index, names=titles)

      return df_mi

  def set_cat_colors(self, cat_colors, axis, cat_index, cat_title=False):
      for inst_ct in cat_colors:
          if cat_title != False:
              cat_name = cat_title + ': ' + inst_ct
          else:
              cat_name = inst_ct

          inst_color = cat_colors[inst_ct]
          self.set_cat_color(axis=axis, cat_index=cat_index, cat_name=cat_name, inst_color=inst_color)

  def set_manual_category(self, col=None, row=None, preferred_cats=None):
    '''
    This method is used to tell Clustergrammer2 that the user wants to define
    a manual category interactively using the dendrogram.
    '''

    self.dat['manual_category'] = {}
    self.dat['manual_category']['col'] = col
    self.dat['manual_category']['row'] = row

    if preferred_cats is not None:
      pref_cats = []
      for inst_row in preferred_cats.index.tolist():
          inst_dict = {}
          inst_dict['name'] = inst_row
          inst_dict['color'] = preferred_cats.loc[inst_row, 'color']
          pref_cats.append(inst_dict)

      if col is not None:
        self.dat['manual_category']['col_cats'] = pref_cats

      if row is not None:
        self.dat['manual_category']['row_cats'] = pref_cats

  def ds_to_original_meta(self, axis):
    clusters = self.meta_ds_col.index.tolist()

    man_cat_title = self.dat['manual_category'][axis]

    for inst_cluster in clusters:

        # get the manual category from the downsampled data
        inst_man_cat = self.meta_ds_col.loc[inst_cluster, man_cat_title]

        # find original labels that are assigned to this cluster
        found_labels = self.meta_col[self.meta_col[self.ds_name] == inst_cluster].index.tolist()

        # update manual category
        self.meta_col.loc[found_labels, man_cat_title] = inst_man_cat

  def get_manual_category(self, tmp):


    for axis in ['col']:

      try:
      ###########################################################

        export_dict = {}

        inst_nodes = self.dat['nodes'][axis]
        # print('get get_manual_category', len(inst_nodes))

        cat_title = self.dat['manual_category'][axis]

        # Category Names
        try:
        ###########################################################

          export_dict[cat_title] = pd.Series(json.loads(self.widget_instance.manual_cat)[axis][cat_title])

          if hasattr(self, 'meta_cat') == True:

            if axis == 'row':
              if self.is_downsampled == False:
                self.meta_row.loc[inst_nodes, cat_title] = export_dict[cat_title]
              else:
                self.meta_ds_row.loc[inst_nodes, cat_title] = export_dict[cat_title]
                self.ds_to_original_meta(axis)

            elif axis == 'col':
              if self.is_downsampled == False:
                # print('> > > > > > > > > > > ')
                # print(export_dict[cat_title])

                # print('.. .. .. .. .. .. .. .. .. .. .. .. ')
                # print(self.meta_col)

                # # for some reason this is breaking when trying to sync metadata
                # # switching to slower row based syncing
                # self.meta_col.loc[inst_nodes, cat_title] = export_dict[cat_title]

                # manually looping through rows in metadata works
                for inst_row in export_dict[cat_title].index.tolist():
                  self.meta_col.loc[inst_row, cat_title] = export_dict[cat_title][inst_row]

                # print('>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>')
                # print(self.meta_col)

              else:
                self.meta_ds_col.loc[inst_nodes, cat_title] = export_dict[cat_title]
                self.ds_to_original_meta(axis)


        except:
          print('unable to load', axis ,' category, please check title')

        # Category Colors
        #######################
        export_dict['cat_colors'] = json.loads(self.widget_instance.manual_cat)['global_cat_colors']

        # # drop title from category colors
        # export_dict['cat_colors'] = {}
        # for inst_cat in ini_new_colors:
        #   if (': ' in inst_cat):
        #     export_dict['cat_colors'][inst_cat.split(': ')[1]] = ini_new_colors[inst_cat]
        #   else:
        #     export_dict['cat_colors'][inst_cat] = ini_new_colors[inst_cat]

        # if hasattr(self, 'meta_cat') == False:
        #   return export_dict

      except:
          # print('failed to parse manual_category')
          pass


  @staticmethod
  def umi_norm(df):
      # umi norm
      barcode_umi_sum = df.sum()
      df_umi = df.div(barcode_umi_sum)
      return df_umi

  @staticmethod
  def make_df_from_cols(cols):
    inst_col = cols[0]

    cat_titles = []
    for inst_info in inst_col[1:]:
        inst_title = inst_info.split(': ')[0]
        cat_titles.append(inst_title)

    clean_cols = []
    for inst_col in cols:
        inst_clean = []
        for inst_info in inst_col:
            if ': ' in inst_info:
                inst_clean.append(inst_info.split(': ')[1])
            else:
                inst_clean.append(inst_info)
        clean_cols.append(tuple(inst_clean))

    df_ini = pd.DataFrame(data=clean_cols).set_index(0)
    mat = df_ini.values
    rows = df_ini.index.tolist()

    df_meta = pd.DataFrame(data=mat, index=rows, columns=cat_titles)

    return df_meta

add_cats(axis, cat_data)

Add categories to rows or columns using cat_data array of objects. Each object in cat_data is a dictionary with one key (category title) and value (rows/column names) that have this category. Categories will be added onto the existing categories and will be added in the order of the objects in the array.

Example cat_data::

[
  {
    "title": "First Category",
    "cats": {
      "true": [
        "ROS1",
        "AAK1"
      ]
    }
  },
  {
    "title": "Second Category",
    "cats": {
      "something": [
        "PDK4"
      ]
    }
  }
]

Source code in src/celldega/clust/__init__.py

def add_cats(self, axis, cat_data):
  '''
  Add categories to rows or columns using cat_data array of objects. Each object in cat_data is a dictionary with one key (category title) and value (rows/column names) that have this category. Categories will be added onto the existing categories and will be added in the order of the objects in the array.

  Example ``cat_data``::


      [
        {
          "title": "First Category",
          "cats": {
            "true": [
              "ROS1",
              "AAK1"
            ]
          }
        },
        {
          "title": "Second Category",
          "cats": {
            "something": [
              "PDK4"
            ]
          }
        }
      ]


  '''
  for inst_data in cat_data:
    categories.add_cats(self, axis, inst_data)

assess_prediction(df_meta, truth, pred)

Generate confusion matrix from y_info

Source code in src/celldega/clust/__init__.py

def assess_prediction(self, df_meta, truth, pred):
    ''' Generate confusion matrix from y_info '''

    y_info = {}
    y_info['true'] = df_meta[truth].values.tolist()
    y_info['pred'] = df_meta[pred].values.tolist()

    sorted_cats = sorted(list(set(y_info['true'] + y_info['pred'])))
    conf_mat = confusion_matrix(y_info['true'], y_info['pred'], labels=sorted_cats)

    # true columns and pred rows
    df_conf = pd.DataFrame(conf_mat, index=sorted_cats, columns=sorted_cats).transpose()

    total_correct = np.trace(df_conf)
    total_pred = df_conf.sum().sum()
    fraction_correct = total_correct/float(total_pred)

    # calculate ser_correct
    correct_list = []
    cat_counts = df_conf.sum(axis=0)
    all_cols = df_conf.columns.tolist()
    for inst_cat in all_cols:
        inst_correct = df_conf.loc[inst_cat, inst_cat] / cat_counts[inst_cat]
        correct_list.append(inst_correct)

    ser_correct = pd.Series(data=correct_list, index=all_cols)

    return df_conf, ser_correct, fraction_correct

clip(lower=None, upper=None)

Trim values at input thresholds using pandas function

Source code in src/celldega/clust/__init__.py

def clip(self, lower=None, upper=None):
  '''
  Trim values at input thresholds using pandas function
  '''
  df = self.export_df()
  df = df.clip(lower=lower, upper=upper)
  self.load_df(df)

cluster(dist_type='cosine', run_clustering=True, dendro=True, views=[], linkage_type='average', sim_mat=False, filter_sim=0.0, calc_cat_pval=False, run_enrichr=None, enrichrgram=None, clust_library='scipy', min_samples=1, min_cluster_size=2)

The main function performs hierarchical clustering, optionally generates filtered views (e.g. row-filtered views), and generates the : visualization_json.

Used to set views equal to ['N_row_sum', 'N_row_var']

Source code in src/celldega/clust/__init__.py

def cluster(self, dist_type='cosine', run_clustering=True,
               dendro=True, views=[],
               linkage_type='average', sim_mat=False, filter_sim=0.0,
               calc_cat_pval=False, run_enrichr=None, enrichrgram=None,
               clust_library='scipy', min_samples=1, min_cluster_size=2):
  '''
  The main function performs hierarchical clustering, optionally generates
  filtered views (e.g. row-filtered views), and generates the :
  ``visualization_json``.

  Used to set views equal to ['N_row_sum', 'N_row_var']
  '''
  initialize_net.viz(self)

  make_clust_fun.make_clust(self, dist_type=dist_type,
                                  run_clustering=run_clustering,
                                  dendro=dendro,
                                  requested_views=views,
                                  linkage_type=linkage_type,
                                  sim_mat=sim_mat,
                                  filter_sim=filter_sim,
                                  calc_cat_pval=calc_cat_pval,
                                  run_enrichr=run_enrichr,
                                  enrichrgram=enrichrgram,
                                  clust_library=clust_library,
                                  min_samples=min_samples,
                                  min_cluster_size=min_cluster_size)

dat_to_df()

Export Pandas DataFrams (will be deprecated).

Source code in src/celldega/clust/__init__.py

def dat_to_df(self):
  '''
  Export Pandas DataFrams (will be deprecated).
  '''
  return data_formats.dat_to_df(self)

df_to_dat(df, define_cat_colors=False)

Load Pandas DataFrame (will be deprecated).

Source code in src/celldega/clust/__init__.py

def df_to_dat(self, df, define_cat_colors=False):
  '''
  Load Pandas DataFrame (will be deprecated).
  '''
  data_formats.df_to_dat(self, df, define_cat_colors)

downsample(df=None, ds_type='kmeans', axis='row', num_samples=100, random_state=1000, ds_name='Downsample', ds_cluster_name='cluster')

Downsample the matrix rows or columns (currently supporting kmeans only). Users can optionally pass in a DataFrame to be downsampled (and this will be incorporated into the network object).

Source code in src/celldega/clust/__init__.py

def downsample(self, df=None, ds_type='kmeans', axis='row', num_samples=100,
               random_state=1000, ds_name='Downsample',
               ds_cluster_name='cluster'):
  '''
  Downsample the matrix rows or columns (currently supporting kmeans only). Users can optionally pass in a DataFrame to be downsampled (and this will be incorporated into the network object).
  '''
  return downsample_fun.main(self, df, ds_type, axis, num_samples, random_state, ds_name, ds_cluster_name)

enrichrgram(lib, axis='row')

Add Enrichr gene enrichment results to your visualization (where your rows are genes). Run enrichrgram before clustering to incldue enrichment results as row categories. Enrichrgram can also be run on the front-end using the Enrichr logo at the top left.

Set lib to the Enrichr library that you want to use for enrichment analysis. Libraries included:

• ChEA_2016
• KEA_2015
• ENCODE_TF_ChIP-seq_2015
• ENCODE_Histone_Modifications_2015
• Disease_Perturbations_from_GEO_up
• Disease_Perturbations_from_GEO_down
• GO_Molecular_Function_2015
• GO_Biological_Process_2015
• GO_Cellular_Component_2015
• Reactome_2016
• KEGG_2016
• MGI_Mammalian_Phenotype_Level_4
• LINCS_L1000_Chem_Pert_up
• LINCS_L1000_Chem_Pert_down

Source code in src/celldega/clust/__init__.py

def enrichrgram(self, lib, axis='row'):
  '''
  Add Enrichr gene enrichment results to your visualization (where your rows
  are genes). Run enrichrgram before clustering to incldue enrichment results
  as row categories. Enrichrgram can also be run on the front-end using the
  Enrichr logo at the top left.

  Set lib to the Enrichr library that you want to use for enrichment analysis.
  Libraries included:

    * ChEA_2016
    * KEA_2015
    * ENCODE_TF_ChIP-seq_2015
    * ENCODE_Histone_Modifications_2015
    * Disease_Perturbations_from_GEO_up
    * Disease_Perturbations_from_GEO_down
    * GO_Molecular_Function_2015
    * GO_Biological_Process_2015
    * GO_Cellular_Component_2015
    * Reactome_2016
    * KEGG_2016
    * MGI_Mammalian_Phenotype_Level_4
    * LINCS_L1000_Chem_Pert_up
    * LINCS_L1000_Chem_Pert_down

  '''

  df = self.export_df()
  df, bar_info = enr_fun.add_enrichr_cats(df, axis, lib)
  self.load_df(df)

  self.dat['enrichrgram_lib'] = lib
  self.dat['row_cat_bars'] = bar_info

export_df()

Export Pandas DataFrame/

Source code in src/celldega/clust/__init__.py

def export_df(self):
  '''
  Export Pandas DataFrame/
  '''
  df = data_formats.dat_to_df(self)

  # drop tuple categories if downsampling
  if self.is_downsampled:
    df.columns = self.dat['nodes']['col']
    df.index = self.dat['nodes']['row']

  return df

export_net_json(net_type='viz', indent='no-indent')

Export dat or viz JSON.

Source code in src/celldega/clust/__init__.py

def export_net_json(self, net_type='viz', indent='no-indent'):
  '''
  Export dat or viz JSON.
  '''
  return export_data.export_net_json(self, net_type, indent)

export_viz_to_widget(which_viz='viz')

Export viz JSON, for use with clustergrammer_widget. Formerly method was named widget.

Source code in src/celldega/clust/__init__.py

def export_viz_to_widget(self, which_viz='viz'):
  '''
  Export viz JSON, for use with clustergrammer_widget. Formerly method was
  named widget.
  '''

  return export_data.export_net_json(self, which_viz, 'no-indent')

filter_N_top(N_top, rank_type='sum', inst_rc=None, axis=None)

Filter the matrix rows or columns based on sum/variance, and only keep the top N.

Source code in src/celldega/clust/__init__.py

def filter_N_top(self, N_top, rank_type='sum', inst_rc=None, axis=None):
  '''
  Filter the matrix rows or columns based on sum/variance, and only keep the top
  N.
  '''

  if axis is None:
    axis = inst_rc
    print('warning inst_rc argument will be deprecated, please use axis')

    if inst_rc is None:
      print('please provide axis argument')

  inst_df = self.dat_to_df()

  inst_df = run_filter.filter_N_top(inst_rc, inst_df, N_top, rank_type)

  self.df_to_dat(inst_df)

filter_cat(axis, cat_index, cat_name)

Filter the matrix based on their category. cat_index is the index of the category, the first category has index=1.

Source code in src/celldega/clust/__init__.py

def filter_cat(self, axis, cat_index, cat_name):
  '''
  Filter the matrix based on their category. cat_index is the index of the category, the first category has index=1.
  '''
  run_filter.filter_cat(self, axis, cat_index, cat_name)

filter_names(axis, names)

Filter the visualization using row/column names. The function takes, axis ('row'/'col') and names, a list of strings.

Source code in src/celldega/clust/__init__.py

def filter_names(self, axis, names):
  '''
  Filter the visualization using row/column names. The function takes, axis ('row'/'col') and names, a list of strings.
  '''
  run_filter.filter_names(self, axis, names)

filter_sum(threshold, take_abs=True, axis=None, inst_rc=None)

Filter a network's rows or columns based on the sum across rows or columns.

Source code in src/celldega/clust/__init__.py

def filter_sum(self, threshold, take_abs=True, axis=None, inst_rc=None):
  '''
  Filter a network's rows or columns based on the sum across rows or columns.
  '''

  if axis is None:
    axis = inst_rc
    print('warning inst_rc argument will be deprecated, please use axis')

    if inst_rc is None:
      print('please provide axis argument')

  inst_df = self.dat_to_df()
  if axis == 'row':
    inst_df = run_filter.df_filter_row_sum(inst_df, threshold, take_abs)
  elif axis == 'col':
    inst_df = run_filter.df_filter_col_sum(inst_df, threshold, take_abs)
  self.df_to_dat(inst_df)

filter_threshold(threshold, num_occur=1, inst_rc=None, axis=None)

Filter the matrix rows or columns based on num_occur values being above a threshold (in absolute value).

Source code in src/celldega/clust/__init__.py

def filter_threshold(self, threshold, num_occur=1, inst_rc=None, axis=None):
  '''
  Filter the matrix rows or columns based on num_occur values being above a
  threshold (in absolute value).
  '''

  if axis is None:
    axis = inst_rc
    print('warning inst_rc argument will be deprecated, please use axis')

    if inst_rc is None:
      print('please provide axis argument')

  inst_df = self.dat_to_df()

  inst_df = run_filter.filter_threshold(inst_df, axis, threshold,
    num_occur)

  self.df_to_dat(inst_df)

generate_signatures(df_data, df_meta, category_name, pval_cutoff=0.05, num_top_dims=False, verbose=True, equal_var=False)

Generate signatures for column categories T-test is run for each category in a one-vs-all manner. The num_top_dims overrides the P-value cutoff.

Source code in src/celldega/clust/__init__.py

def generate_signatures(self, df_data, df_meta, category_name, pval_cutoff=0.05,
                        num_top_dims=False, verbose=True, equal_var=False):

    '''
    Generate signatures for column categories
    T-test is run for each category in a one-vs-all manner. The num_top_dims overrides
    the P-value cutoff.
    '''

    df_t = df_data.transpose()

    # remove columns (dimensions) with constant values
    orig_num_cols = df_t.shape[1]
    df_t = df_t.loc[:, (df_t != df_t.iloc[0]).any()]
    if df_t.shape[1] < orig_num_cols:
      print('dropped columns with constant values')

    # make tuple rows
    df_t.index = [(x, df_meta.loc[x, category_name]) for x in df_t.index.tolist()]
    category_level = 1

    df = self.row_tuple_to_multiindex(df_t)

    cell_types = sorted(list(set(df.index.get_level_values(category_level).tolist())))

    keep_genes = []
    keep_genes_dict = {}
    gene_pval_dict = {}
    all_fold_info = {}

    for inst_ct in cell_types:

        inst_ct_mat = df.xs(key=inst_ct, level=category_level)
        inst_other_mat = df.drop(inst_ct, level=category_level)

        # save mean values and fold change
        fold_info = {}
        fold_info['cluster_mean'] = inst_ct_mat.mean()
        fold_info['other_mean'] = inst_other_mat.mean()
        fold_info['log2_fold'] = fold_info['cluster_mean']/fold_info['other_mean']
        fold_info['log2_fold'] = fold_info['log2_fold'].apply(np.log2)
        all_fold_info[inst_ct] = fold_info

        inst_stats, inst_pvals = ttest_ind(inst_ct_mat, inst_other_mat, axis=0, equal_var=equal_var)

        ser_pval = pd.Series(data=inst_pvals, index=df.columns.tolist()).sort_values()

        if num_top_dims == False:
            ser_pval_keep = ser_pval[ser_pval < pval_cutoff]
        else:
            ser_pval_keep = ser_pval[:num_top_dims]

        gene_pval_dict[inst_ct] = ser_pval_keep

        inst_keep = ser_pval_keep.index.tolist()
        keep_genes.extend(inst_keep)
        keep_genes_dict[inst_ct] = inst_keep

    keep_genes = sorted(list(set(keep_genes)))

    df_gbm = df.groupby(level=category_level).mean().transpose()
    cols = df_gbm.columns.tolist()

    df_sig = df_gbm.loc[keep_genes]

    if len(keep_genes) == 0 and verbose:
        print('found no informative dimensions')

    df_gene_pval = pd.concat(gene_pval_dict, axis=1, sort=False)

    df_diff = {}
    for inst_col in df_gene_pval:

        inst_pvals = df_gene_pval[inst_col]

        # nans represent dimensions that did not meet pval or top threshold
        inst_pvals = inst_pvals.dropna().sort_values()

        inst_genes = inst_pvals.index.tolist()

        # prevent failure if no cells have this category
        if inst_pvals.shape[0] > 0:
            rej, pval_corr = smm.multipletests(inst_pvals, 0.05, method='fdr_bh')[:2]

            ser_pval = pd.Series(data=inst_pvals, index=inst_genes, name='P-values')
            ser_bh_pval = pd.Series(data=pval_corr, index=inst_genes, name='BH P-values')
            ser_log2_fold = all_fold_info[inst_col]['log2_fold'].loc[inst_genes]
            ser_log2_fold.name = 'Log2 Fold Change'
            ser_cluster_mean = all_fold_info[inst_col]['cluster_mean'].loc[inst_genes]
            ser_cluster_mean.name = 'Cluster Mean'
            ser_other_mean = all_fold_info[inst_col]['other_mean'].loc[inst_genes]
            ser_other_mean.name = 'All Other Mean'
            inst_df = pd.concat([ser_pval, ser_bh_pval, ser_log2_fold, ser_cluster_mean, ser_other_mean], axis=1)

            df_diff[inst_col] = inst_df

    return df_sig, df_diff

load_data_file_to_net(filename)

Load Clustergrammer's dat format (saved as JSON).

Source code in src/celldega/clust/__init__.py

def load_data_file_to_net(self, filename):
  '''
  Load Clustergrammer's dat format (saved as JSON).
  '''
  inst_dat = self.load_json_to_dict(filename)
  load_data.load_data_to_net(self, inst_dat)

load_df(df_ini, meta_col=None, meta_row=None, col_cats=None, row_cats=None, is_downsampled=False, meta_ds_row=None, meta_ds_col=None)

Load Pandas DataFrame.

Source code in src/celldega/clust/__init__.py

def load_df(self, df_ini, meta_col=None, meta_row=None, col_cats=None,
            row_cats=None, is_downsampled=False, meta_ds_row=None,
            meta_ds_col=None):
  '''
  Load Pandas DataFrame.
  '''
  self.reset()

  # load dataframe
  df = deepcopy(df_ini)


  if is_downsampled:
    if meta_ds_col is not None:
      self.meta_ds_col = meta_ds_col
    if meta_ds_row is not None:
      self.meta_ds_row = meta_ds_row

  # define downsampled status
  self.is_downsampled = is_downsampled
  # print('load_df: is_downsampled', is_downsampled)

  if hasattr(self, 'meta_col') == False and hasattr(self, 'meta_row') == False:
    self.meta_cat = False

  # load metadata
  if isinstance(meta_col, pd.DataFrame):
    self.meta_col = meta_col

    if col_cats is None:
      self.col_cats = meta_col.columns.tolist()
    else:
      self.col_cats = col_cats

    self.meta_cat = True

  if isinstance(meta_row, pd.DataFrame):
    self.meta_row = meta_row

    if row_cats is None:
      self.row_cats = meta_row.columns.tolist()
    else:
      self.row_cats = row_cats

    self.meta_cat = True

  data_formats.df_to_dat(self, df, define_cat_colors=True)

load_file(filename)

Load TSV file.

Source code in src/celldega/clust/__init__.py

def load_file(self, filename):
  '''
  Load TSV file.
  '''
  load_data.load_file(self, filename)

load_file_as_string(file_string, filename='')

Load file as a string.

Source code in src/celldega/clust/__init__.py

def load_file_as_string(self, file_string, filename=''):
  '''
  Load file as a string.
  '''
  load_data.load_file_as_string(self, file_string, filename=filename)

load_gene_exp_to_df(inst_path) staticmethod

Loads gene expression data from 10x in sparse matrix format and returns a Pandas dataframe

Source code in src/celldega/clust/__init__.py

@staticmethod
def load_gene_exp_to_df(inst_path):
  '''
  Loads gene expression data from 10x in sparse matrix format and returns a
  Pandas dataframe
  '''

  import pandas as pd
  from scipy import io
  from scipy import sparse
  from ast import literal_eval as make_tuple

  # matrix
  Matrix = io.mmread( inst_path + 'matrix.mtx')
  mat = Matrix.todense()

  # genes
  filename = inst_path + 'genes.tsv'
  f = open(filename, 'r')
  lines = f.readlines()
  f.close()

  # # add unique id to all genes
  # genes = []
  # unique_id = 0
  # for inst_line in lines:
  #     inst_line = inst_line.strip().split()

  #     if len(inst_line) > 1:
  #       inst_gene = inst_line[1]
  #     else:
  #       inst_gene = inst_line[0]

  #     genes.append(inst_gene + '_' + str(unique_id))
  #     unique_id = unique_id + 1

  # add unique id only to duplicate genes
  ini_genes = []
  for inst_line in lines:
      inst_line = inst_line.strip().split()
      if len(inst_line) > 1:
        inst_gene = inst_line[1]
      else:
        inst_gene = inst_line[0]
      ini_genes.append(inst_gene)

  gene_name_count = pd.Series(ini_genes).value_counts()
  duplicate_genes = gene_name_count[gene_name_count > 1].index.tolist()

  dup_index = {}
  genes = []
  for inst_row in ini_genes:

    # add index to non-unique genes
    if inst_row in duplicate_genes:

      # calc_non-unque index
      if inst_row not in dup_index:
        dup_index[inst_row] = 1
      else:
        dup_index[inst_row] = dup_index[inst_row] + 1

      new_row = inst_row + '_' + str(dup_index[inst_row])

    else:
      new_row = inst_row

    genes.append(new_row)

  # barcodes
  filename = inst_path + 'barcodes.tsv'
  f = open(filename, 'r')
  lines = f.readlines()
  f.close()

  cell_barcodes = []
  for inst_bc in lines:
      inst_bc = inst_bc.strip().split('\t')

      # remove dash from barcodes if necessary
      if '-' in inst_bc[0]:
        inst_bc[0] = inst_bc[0].split('-')[0]

      cell_barcodes.append(inst_bc[0])

  # parse tuples if necessary
  try:
      cell_barcodes = [make_tuple(x) for x in cell_barcodes]
  except:
      pass

  try:
      genes = [make_tuple(x) for x in genes]
  except:
      pass

  # make dataframe
  df = pd.DataFrame(mat, index=genes, columns=cell_barcodes)

  return df

load_stdin()

Load stdin TSV-formatted string.

Source code in src/celldega/clust/__init__.py

def load_stdin(self):
  '''
  Load stdin TSV-formatted string.
  '''
  load_data.load_stdin(self)

load_tsv_to_net(file_buffer, filename=None)

This will load a TSV matrix file buffer; this is exposed so that it will be possible to load data without having to read from a file.

Source code in src/celldega/clust/__init__.py

def load_tsv_to_net(self, file_buffer, filename=None):
  '''
  This will load a TSV matrix file buffer; this is exposed so that it will
  be possible to load data without having to read from a file.
  '''
  load_data.load_tsv_to_net(self, file_buffer, filename)

load_vect_post_to_net(vect_post)

Load data in the vector format JSON.

Source code in src/celldega/clust/__init__.py

def load_vect_post_to_net(self, vect_post):
  '''
  Load data in the vector format JSON.
  '''
  load_vect_post.main(self, vect_post)

normalize(df=None, norm_type='zscore', axis='row', z_clip=None)

Normalize the matrix rows or columns using Z-score (zscore) or Quantile Normalization (qn). Users can optionally pass in a DataFrame to be normalized (and this will be incorporated into the Network object).

Source code in src/celldega/clust/__init__.py

def normalize(self, df=None, norm_type='zscore', axis='row', z_clip=None):
  '''
  Normalize the matrix rows or columns using Z-score (zscore) or Quantile Normalization (qn). Users can optionally pass in a DataFrame to be normalized (and this will be incorporated into the Network object).
  '''
  normalize_fun.run_norm(self, df, norm_type, axis, z_clip)

old_confusion_matrix_and_correct_series(y_info)

Generate confusion matrix from y_info

Source code in src/celldega/clust/__init__.py

def old_confusion_matrix_and_correct_series(self, y_info):
    ''' Generate confusion matrix from y_info '''


    a = deepcopy(y_info['true'])
    true_count = dict((i, a.count(i)) for i in set(a))

    a = deepcopy(y_info['pred'])
    pred_count = dict((i, a.count(i)) for i in set(a))

    sorted_cats = sorted(list(set(y_info['true'] + y_info['pred'])))
    conf_mat = confusion_matrix(y_info['true'], y_info['pred'], sorted_cats)
    df_conf = pd.DataFrame(conf_mat, index=sorted_cats, columns=sorted_cats)

    total_correct = np.trace(df_conf)
    total_pred = df_conf.sum().sum()
    fraction_correct = total_correct/float(total_pred)

    # calculate ser_correct
    correct_list = []
    cat_counts = df_conf.sum(axis=1)
    all_cols = df_conf.columns.tolist()
    for inst_cat in all_cols:
        inst_correct = df_conf[inst_cat].loc[inst_cat] / cat_counts[inst_cat]
        correct_list.append(inst_correct)

    ser_correct = pd.Series(data=correct_list, index=all_cols)

    populations = {}
    populations['true'] = true_count
    populations['pred'] = pred_count

    return df_conf, populations, ser_correct, fraction_correct

old_generate_signatures(df_ini, category_level, pval_cutoff=0.05, num_top_dims=False, verbose=True, equal_var=False)

Generate signatures for column categories

Source code in src/celldega/clust/__init__.py

def old_generate_signatures(self, df_ini, category_level, pval_cutoff=0.05,
                        num_top_dims=False, verbose=True, equal_var=False):

    ''' Generate signatures for column categories '''

    # change this to use category name and set caetgory level to 1 always
    ###################################

    df_t = df_ini.transpose()

    # remove columns with constant values
    df_t = df_t.loc[:, (df_t != df_t.iloc[0]).any()]

    df = self.row_tuple_to_multiindex(df_t)

    cell_types = sorted(list(set(df.index.get_level_values(category_level).tolist())))

    keep_genes = []
    keep_genes_dict = {}
    gene_pval_dict = {}
    all_fold_info = {}

    for inst_ct in cell_types:

        inst_ct_mat = df.xs(key=inst_ct, level=category_level)
        inst_other_mat = df.drop(inst_ct, level=category_level)

        # save mean values and fold change
        fold_info = {}
        fold_info['cluster_mean'] = inst_ct_mat.mean()
        fold_info['other_mean'] = inst_other_mat.mean()
        fold_info['log2_fold'] = fold_info['cluster_mean']/fold_info['other_mean']
        fold_info['log2_fold'] = fold_info['log2_fold'].apply(np.log2)
        all_fold_info[inst_ct] = fold_info

        inst_stats, inst_pvals = ttest_ind(inst_ct_mat, inst_other_mat, axis=0, equal_var=equal_var)

        ser_pval = pd.Series(data=inst_pvals, index=df.columns.tolist()).sort_values()

        if num_top_dims == False:
            ser_pval_keep = ser_pval[ser_pval < pval_cutoff]
        else:
            ser_pval_keep = ser_pval[:num_top_dims]

        gene_pval_dict[inst_ct] = ser_pval_keep

        inst_keep = ser_pval_keep.index.tolist()
        keep_genes.extend(inst_keep)
        keep_genes_dict[inst_ct] = inst_keep

    keep_genes = sorted(list(set(keep_genes)))

    df_gbm = df.groupby(level=category_level).mean().transpose()
    cols = df_gbm.columns.tolist()
    new_cols = []
    for inst_col in cols:
        new_col = (inst_col, category_level + ': ' + inst_col)
        new_cols.append(new_col)
    df_gbm.columns = new_cols

    df_sig = df_gbm.loc[keep_genes]

    if len(keep_genes) == 0 and verbose:
        print('found no informative dimensions')

    df_gene_pval = pd.concat(gene_pval_dict, axis=1, sort=False)

    return df_sig, df_gene_pval, all_fold_info

old_predict_cats_from_sigs(df_data_ini, df_sig_ini, dist_type='cosine', predict_level='Predict Category', truth_level=1, unknown_thresh=-1)

Predict category using signature

Source code in src/celldega/clust/__init__.py

def old_predict_cats_from_sigs(self, df_data_ini, df_sig_ini, dist_type='cosine', predict_level='Predict Category',
                           truth_level=1, unknown_thresh=-1):
    ''' Predict category using signature '''

    keep_rows = df_sig_ini.index.tolist()
    data_rows = df_data_ini.index.tolist()

    common_rows = list(set(data_rows).intersection(keep_rows))

    df_data = deepcopy(df_data_ini.loc[common_rows])
    df_sig = deepcopy(df_sig_ini.loc[common_rows])

    # calculate sim_mat of df_data and df_sig
    cell_types = df_sig.columns.tolist()
    barcodes = df_data.columns.tolist()
    sim_mat = 1 - pairwise_distances(df_sig.transpose(), df_data.transpose(), metric=dist_type)
    df_sim = pd.DataFrame(data=sim_mat, index=cell_types, columns=barcodes).transpose()

    # get the top column value (most similar signature)
    df_sim_top = df_sim.idxmax(axis=1)

    # get the maximum similarity of a cell to a cell type definition
    max_sim = df_sim.max(axis=1)

    unknown_cells = max_sim[max_sim < unknown_thresh].index.tolist()

    # assign unknown cells (need category of same name)
    df_sim_top[unknown_cells] = 'Unknown'

    # add predicted category name to top list
    top_list = df_sim_top.values
    top_list = [ predict_level + ': ' + x[0] if type(x) is tuple else predict_level + ': ' + x  for x in top_list]

    # add cell type category to input data
    df_cat = deepcopy(df_data)
    cols = df_cat.columns.tolist()
    new_cols = []

    # check whether the columns have the true category available
    has_truth = False
    if type(cols[0]) is tuple:
        has_truth = True

    if has_truth:
        new_cols = [tuple(list(a) + [b]) for a,b in zip(cols, top_list)]
    else:
        new_cols = [tuple([a] + [b]) for a,b in zip(cols, top_list)]

    # transfer new categories
    df_cat.columns = new_cols

    # keep track of true and predicted labels
    y_info = {}
    y_info['true'] = []
    y_info['pred'] = []

    if has_truth:
        y_info['true'] = [x[truth_level].split(': ')[1] for x in cols]
        y_info['pred'] = [x.split(': ')[1] for x in top_list]

    return df_cat, df_sim.transpose(), y_info

predict_cats_from_sigs(df_data_ini, df_meta, df_sig_ini, predict='Predicted Category', dist_type='cosine', unknown_thresh=-1)

Predict category using signature

Source code in src/celldega/clust/__init__.py

def predict_cats_from_sigs(self, df_data_ini, df_meta, df_sig_ini,
                           predict='Predicted Category', dist_type='cosine',
                           unknown_thresh=-1):
    ''' Predict category using signature '''

    keep_rows = df_sig_ini.index.tolist()
    data_rows = df_data_ini.index.tolist()

    common_rows = list(set(data_rows).intersection(keep_rows))

    df_data = deepcopy(df_data_ini.loc[common_rows])
    df_sig = deepcopy(df_sig_ini.loc[common_rows])

    # calculate sim_mat of df_data and df_sig
    cell_types = df_sig.columns.tolist()
    barcodes = df_data.columns.tolist()
    sim_mat = 1 - pairwise_distances(df_sig.transpose(), df_data.transpose(), metric=dist_type)
    df_sim = pd.DataFrame(data=sim_mat, index=cell_types, columns=barcodes).transpose()

    # get the top column value (most similar signature)
    df_sim_top = df_sim.idxmax(axis=1)

    # get the maximum similarity of a cell to a cell type definition
    max_sim = df_sim.max(axis=1)

    unknown_cells = max_sim[max_sim < unknown_thresh].index.tolist()

    # assign unknown cells (need category of same name)
    df_sim_top[unknown_cells] = 'Unknown'

    # add predicted category name to top list
    # top_list = df_sim_top.values
    df_sim = df_sim.transpose()

    df_meta[predict] = df_sim_top

    return df_sim, df_meta

random_sample(num_samples, df=None, replace=False, weights=None, random_state=100, axis='row')

Return random sample of matrix.

Source code in src/celldega/clust/__init__.py

def random_sample(self, num_samples, df=None, replace=False, weights=None, random_state=100, axis='row'):
  '''
  Return random sample of matrix.
  '''

  if df is None:
    df = self.dat_to_df()

  if axis == 'row':
    axis = 0
  if axis == 'col':
    axis = 1

  df = self.export_df()
  df = df.sample(n=num_samples, replace=replace, weights=weights, random_state=random_state,  axis=axis)

  self.load_df(df)

reset()

This re-initializes the Network object.

Source code in src/celldega/clust/__init__.py

def reset(self):
  '''
  This re-initializes the Network object.
  '''
  initialize_net.main(self)

set_manual_category(col=None, row=None, preferred_cats=None)

This method is used to tell Clustergrammer2 that the user wants to define a manual category interactively using the dendrogram.

Source code in src/celldega/clust/__init__.py

def set_manual_category(self, col=None, row=None, preferred_cats=None):
  '''
  This method is used to tell Clustergrammer2 that the user wants to define
  a manual category interactively using the dendrogram.
  '''

  self.dat['manual_category'] = {}
  self.dat['manual_category']['col'] = col
  self.dat['manual_category']['row'] = row

  if preferred_cats is not None:
    pref_cats = []
    for inst_row in preferred_cats.index.tolist():
        inst_dict = {}
        inst_dict['name'] = inst_row
        inst_dict['color'] = preferred_cats.loc[inst_row, 'color']
        pref_cats.append(inst_dict)

    if col is not None:
      self.dat['manual_category']['col_cats'] = pref_cats

    if row is not None:
      self.dat['manual_category']['row_cats'] = pref_cats

sim_same_and_diff_category_samples(df, cat_index=1, dist_type='cosine', equal_var=False, plot_roc=True, precalc_dist=False, calc_roc=True)

Calculate the similarity of samples from the same and different categories. The cat_index gives the index of the category, where 1 in the first category

Source code in src/celldega/clust/__init__.py

def sim_same_and_diff_category_samples(self, df, cat_index=1, dist_type='cosine',
                                       equal_var=False, plot_roc=True,
                                       precalc_dist=False, calc_roc=True):
    '''
    Calculate the similarity of samples from the same and different categories. The
    cat_index gives the index of the category, where 1 in the first category
    '''

    cols = df.columns.tolist()

    if type(precalc_dist) == bool:
        # compute distnace between rows (transpose to get cols as rows)
        dist_arr = 1 - pdist(df.transpose(), metric=dist_type)
    else:
        dist_arr = precalc_dist

    # generate sample names with categories
    sample_combos = list(combinations(range(df.shape[1]),2))

    sample_names = [str(ind) + '_same' if cols[x[0]][cat_index] == cols[x[1]][cat_index] else str(ind) + '_different' for ind, x in enumerate(sample_combos)]

    ser_dist = pd.Series(data=dist_arr, index=sample_names)

    # find same-cat sample comparisons
    same_cat = [x for x in sample_names if x.split('_')[1] == 'same']

    # find diff-cat sample comparisons
    diff_cat = [x for x in sample_names if x.split('_')[1] == 'different']

    # make series of same and diff category sample comparisons
    ser_same = ser_dist[same_cat]
    ser_same.name = 'Same Category'
    ser_diff = ser_dist[diff_cat]
    ser_diff.name = 'Different Category'

    sim_dict = {}
    roc_data = {}
    sim_data = {}

    sim_dict['same'] = ser_same
    sim_dict['diff'] = ser_diff

    pval_dict = {}
    ttest_stat, pval_dict['ttest'] = ttest_ind(ser_diff, ser_same, equal_var=equal_var)

    ttest_stat, pval_dict['mannwhitney'] = mannwhitneyu(ser_diff, ser_same)

    if calc_roc:
        # calc AUC
        true_index = list(np.ones(sim_dict['same'].shape[0]))
        false_index = list(np.zeros(sim_dict['diff'].shape[0]))
        y_true = true_index + false_index

        true_val = list(sim_dict['same'].values)
        false_val = list(sim_dict['diff'].values)
        y_score = true_val + false_val

        fpr, tpr, thresholds = roc_curve(y_true, y_score)

        inst_auc = auc(fpr, tpr)

        if plot_roc:
            plt.figure()
            plt.plot(fpr, tpr)
            plt.plot([0, 1], [0, 1], color='navy', linestyle='--')
            plt.figure(figsize=(10,10))

            print('AUC', inst_auc)

        roc_data['true'] = y_true
        roc_data['score'] = y_score
        roc_data['fpr'] = fpr
        roc_data['tpr'] = tpr
        roc_data['thresholds'] = thresholds
        roc_data['auc'] = inst_auc

    sim_data['sim_dict'] = sim_dict
    sim_data['pval_dict'] = pval_dict
    sim_data['roc_data'] = roc_data

    return sim_data

swap_nan_for_zero()

Swaps all NaN (numpy NaN) instances for zero.

Source code in src/celldega/clust/__init__.py

def swap_nan_for_zero(self):
  '''
  Swaps all NaN (numpy NaN) instances for zero.
  '''
  self.dat['mat'][np.isnan(self.dat['mat'])] = 0

widget(which_viz='viz', link_net=None, link_net_js=None, clust_library='scipy', min_samples=1, min_cluster_size=2)

Generate a widget visualization using the widget. The export_viz_to_widget method passes the visualization JSON to the instantiated widget, which is returned and visualized on the front-end.

Source code in src/celldega/clust/__init__.py

def widget(self, which_viz='viz', link_net=None, link_net_js=None, clust_library='scipy',
  min_samples=1, min_cluster_size=2):
  '''
  Generate a widget visualization using the widget. The export_viz_to_widget
  method passes the visualization JSON to the instantiated widget, which is
  returned and visualized on the front-end.
  '''
  # run clustering if necessary
  if len(self.viz['row_nodes']) == 0:
    self.cluster(clust_library=clust_library, min_samples=min_samples,
                 min_cluster_size=min_cluster_size)

    # add manual_category to viz json
    if 'manual_category' in self.dat:
      self.viz['manual_category'] = self.dat['manual_category']

    # add pre-z-score data to viz
    if 'pre_zscore' in self.dat:
      self.viz['pre_zscore'] = self.dat['pre_zscore']

  self.widget_instance = self.widget_class(network = self.export_viz_to_widget(which_viz))

  # initialize manual category
  if 'manual_category' in self.dat:
    manual_cat = {}
    axis = 'col'
    manual_cat[axis] = {}
    manual_cat[axis]['col_cat_colors'] = self.viz['cat_colors'][axis]['cat-0']

    man_cat_name = self.dat['manual_category'][axis]
    if axis == 'col':
      manual_cat[axis][man_cat_name] = self.meta_col[man_cat_name].to_dict()

    self.widget_instance.manual_cat = json.dumps(manual_cat)

    self.widget_instance.observe(self.get_manual_category, names='manual_cat')

  # add link (python)
  if link_net is not None:
    inst_link = widgets.link(
                              (self.widget_instance, 'manual_cat'),
                              (link_net.widget_instance, 'manual_cat')
                             )
    self.widget_instance.link = inst_link

  # add jslink (JavaScript)
  if link_net_js is not None:

    inst_link = widgets.jslink(
                              (self.widget_instance, 'manual_cat'),
                              (link_net_js.widget_instance, 'manual_cat')
                             )
    self.widget_instance.link = inst_link

  return self.widget_instance

widget_df()

Export a DataFrame from the front-end visualization. For instance, a user can filter to show only a single cluster using the dendrogram and then get a dataframe of this cluster using the widget_df method.

Source code in src/celldega/clust/__init__.py

def widget_df(self):
  '''
  Export a DataFrame from the front-end visualization. For instance, a user
  can filter to show only a single cluster using the dendrogram and then
  get a dataframe of this cluster using the widget_df method.
  '''

  if hasattr(self, 'widget_instance') == True:

    if self.widget_instance.mat_string != '':

      tmp_net = deepcopy(Network())

      df_string = self.widget_instance.mat_string

      tmp_net.load_file_as_string(df_string)

      df = tmp_net.export_df()

      return df

    else:
      return self.export_df()

  else:
    if hasattr(self, 'widget_class') == True:
      print('Please make the widget before exporting the widget DataFrame.')
      print('Do this using the widget method: net.widget()')

    else:
      print('Can not make widget because Network has no attribute widget_class')
      print('Please instantiate Network with clustergrammer_widget using: Network(clustergrammer_widget)')

write_json_to_file(net_type, filename, indent='no-indent')

Save dat or viz as a JSON to file.

Source code in src/celldega/clust/__init__.py

def write_json_to_file(self, net_type, filename, indent='no-indent'):
  '''
  Save dat or viz as a JSON to file.
  '''
  export_data.write_json_to_file(self, net_type, filename, indent)

write_matrix_to_tsv(filename=None, df=None)

Export data-matrix to file.

Source code in src/celldega/clust/__init__.py

def write_matrix_to_tsv(self, filename=None, df=None):
  '''
  Export data-matrix to file.
  '''
  return export_data.write_matrix_to_tsv(self, filename, df)

batch_transform_geometries(geometries, transformation_matrix, scale)

Batch transform geometries using numpy for optimized performance.

Source code in src/celldega/pre/boundary_tile.py

def batch_transform_geometries(geometries, transformation_matrix, scale):
    """
    Batch transform geometries using numpy for optimized performance.
    """
    # Extract affine transformation parameters into a 3x3 matrix for numpy
    affine_matrix = np.array(
        [
            [
                transformation_matrix[0, 0],
                transformation_matrix[0, 1],
                transformation_matrix[0, 2],
            ],
            [
                transformation_matrix[1, 0],
                transformation_matrix[1, 1],
                transformation_matrix[1, 2],
            ],
            [0, 0, 1],
        ]
    )

    transformed_geometries = []

    for geom in geometries:

        if isinstance(geom, Point):
            # Transform a single Point geometry
            point_coords = np.array([[geom.x, geom.y]])
            transformed_coords = (
                numpy_affine_transform(point_coords, affine_matrix) / scale
            )
            transformed_geometries.append(Point(transformed_coords[0]))

        elif isinstance(geom, Polygon):
            # Transform a Polygon geometry
            exterior_coords = np.array(geom.exterior.coords)

            # Apply the affine transformation and scale
            transformed_coords = (
                numpy_affine_transform(exterior_coords, affine_matrix) / scale
            )

            # Append the result to the transformed_geometries list
            transformed_geometries.append([transformed_coords.tolist()])

        elif isinstance(geom, MultiPolygon):
            geom = next(geom.geoms)  # Use the first geometry

            # Transform the exterior of the polygon
            exterior_coords = np.array(geom.exterior.coords)

            # Apply the affine transformation and scale
            transformed_coords = (
                numpy_affine_transform(exterior_coords, affine_matrix) / scale
            )

            # Append the result to the transformed_geometries list
            transformed_geometries.append([transformed_coords.tolist()])

    return transformed_geometries

calc_meta_gene_data(cbg)

Calculate gene metadata from the cell-by-gene matrix

Parameters:

Name	Type	Description	Default
cbg	DataFrame	A sparse DataFrame with genes as columns and barcodes as rows.	required

Returns:

Type	Description
	pandas.DataFrame: A DataFrame with gene metadata including mean, standard deviation, maximum expression, and proportion of non-zero expression.

Source code in src/celldega/pre/landscape.py

def calc_meta_gene_data(cbg):
    """
    Calculate gene metadata from the cell-by-gene matrix

    Args:
        cbg (pandas.DataFrame): A sparse DataFrame with genes as columns and barcodes as rows.

    Returns:
        pandas.DataFrame: A DataFrame with gene metadata including mean, standard deviation,
            maximum expression, and proportion of non-zero expression.
    """

    # Helper function to convert to dense if sparse
    def convert_to_dense(series):
        """
        Convert a pandas Series to dense format if it's sparse.

        Parameters
        ----------
        series : pandas.Series

        Returns
        -------
        pandas.Series
            Dense Series if input was sparse; original Series otherwise.
        """
        if pd.api.types.is_sparse(series):
            return series.sparse.to_dense()
        return series

    # Ensure cbg is a DataFrame
    if not isinstance(cbg, pd.DataFrame):
        raise TypeError("cbg must be a pandas DataFrame")

    # Determine if cbg is sparse
    is_sparse = pd.api.types.is_sparse(cbg)

    if is_sparse:
        # Ensure cbg has SparseDtype with float and fill_value=0
        cbg = cbg.astype(pd.SparseDtype("float", fill_value=0))
        print("cbg is a sparse DataFrame. Proceeding with sparse operations.")
    else:
        print("cbg is a dense DataFrame. Proceeding with dense operations.")

    # Calculate mean expression across tiles
    print("Calculating mean expression")
    mean_expression = cbg.mean(axis=0)

    # Calculate variance as the average of the squared deviations
    print("Calculating variance")
    num_tiles = cbg.shape[1]
    variance = cbg.apply(
        lambda x: ((x - mean_expression[x.name]) ** 2).sum() / num_tiles, axis=0
    )
    std_deviation = np.sqrt(variance)

    # Calculate maximum expression
    max_expression = cbg.max(axis=0)

    # Calculate proportion of tiles with non-zero expression
    proportion_nonzero = (cbg != 0).sum(axis=0) / len(cbg)

    # Create a DataFrame to hold all these metrics

    meta_gene = pd.DataFrame(
        {
            "mean": mean_expression.sparse.to_dense() if isinstance(mean_expression.dtype, pd.SparseDtype) else mean_expression,
            "std": std_deviation,
            "max": max_expression.sparse.to_dense() if isinstance(max_expression.dtype, pd.SparseDtype) else max_expression,
            "non-zero": proportion_nonzero.sparse.to_dense() if isinstance(proportion_nonzero.dtype, pd.SparseDtype) else proportion_nonzero,
        }
    )

    meta_gene_clean = pd.DataFrame(
        meta_gene.values, index=meta_gene.index.tolist(), columns=meta_gene.columns
    )

    return meta_gene_clean

cluster_gene_expression(technology, path_landscape_files, cbg, data_dir=None, segmentation_approach='default')

Calculates cluster-specific gene expression signatures for Xenium data.

Parameters:

Name	Type	Description	Default
technology	str	The technology used (e.g., "Xenium" or "MERSCOPE"). Currently, only "Xenium" is supported.	required
data_dir	str	Path to the directory containing the Xenium data.	None
path_landscape_files	str	Path to the directory where the gene expression signature file will be saved.	required
cbg	DataFrame	A cell-by-gene matrix where rows represent cells and columns represent genes. The index of the DataFrame should match the cell IDs in the Xenium metadata.	required

Raises:

Type	Description
ValueError	If the specified technology is not supported.
FileNotFoundError	If the required input files are not found.

Source code in src/celldega/pre/__init__.py

def cluster_gene_expression(technology, path_landscape_files, cbg, data_dir=None, segmentation_approach='default'):
    """
    Calculates cluster-specific gene expression signatures for Xenium data.

    Args:
        technology (str): The technology used (e.g., "Xenium" or "MERSCOPE"). Currently, only "Xenium" is supported.
        data_dir (str): Path to the directory containing the Xenium data.
        path_landscape_files (str): Path to the directory where the gene expression signature file will be saved.
        cbg (pd.DataFrame): A cell-by-gene matrix where rows represent cells and columns represent genes.
                            The index of the DataFrame should match the cell IDs in the Xenium metadata.

    Raises:
        ValueError: If the specified technology is not supported.
        FileNotFoundError: If the required input files are not found.
    """

    print("\n========Create cluster gene expression (df_sig)========")
    if technology != "Xenium" and technology != 'custom':
        raise ValueError(
            f"Unsupported technology: {technology}. Currently, only 'Xenium' and 'Custom' is supported."
        )

    if technology == "Xenium":
        cells_csv_path = f'{data_dir}/cells.csv.gz'
        clusters_csv_path = (
            f'{data_dir}/analysis/clustering/gene_expression_graphclust/clusters.csv'
        )

        # Load the cell metadata
        usecols = ['cell_id', 'x_centroid', 'y_centroid']
        meta_cell = pd.read_csv(cells_csv_path, index_col=0, usecols=usecols)
        meta_cell.columns = ['center_x', 'center_y']

        # Load the clustering data
        df_meta = pd.read_csv(clusters_csv_path, index_col=0)
        df_meta['Cluster'] = df_meta['Cluster'].astype('string')
        df_meta.columns = ['cluster']

        # Add cluster information to the cell metadata
        meta_cell['cluster'] = df_meta['cluster']
        clusters = meta_cell['cluster'].unique().tolist()

        # Calculate cluster-specific gene expression signatures
        list_ser = []
        for inst_cat in meta_cell['cluster'].unique().tolist():
            if inst_cat is not None:
                inst_cells = meta_cell[meta_cell['cluster'] == inst_cat].index.tolist()
                inst_ser = cbg.loc[inst_cells].sum() / len(inst_cells)
                inst_ser.name = inst_cat
                list_ser.append(inst_ser)

    elif technology == "custom":

        df_cluster = pd.read_parquet(os.path.join(path_landscape_files, f"cell_clusters_{segmentation_approach}", "cluster.parquet"))
        clusters = df_cluster['cluster'].unique().tolist()

        list_ser = []
        for inst_cat in df_cluster['cluster'].unique():
            if inst_cat is not None:
                inst_cells = df_cluster[df_cluster['cluster'] == inst_cat].index.tolist()

                if set(inst_cells) & set(cbg.index):
                    common_cells = list(set(inst_cells) & set(cbg.index))
                    inst_ser = cbg.loc[common_cells].sum()/len(common_cells)

                else:
                    genes = cbg.columns
                    inst_ser = pd.Series(0.0, index=genes)

                inst_ser.name = inst_cat
                list_ser.append(inst_ser)

    # Combine the signatures into a DataFrame
    df_sig = pd.concat(list_ser, axis=1)

    # Handle potential multiindex issues
    df_sig.columns = df_sig.columns.tolist()
    df_sig.index = df_sig.index.tolist()

    # Filter out unwanted genes
    keep_genes = df_sig.index.tolist()
    keep_genes = [x for x in keep_genes if 'Unassigned' not in x]
    keep_genes = [x for x in keep_genes if 'NegControl' not in x]
    keep_genes = [x for x in keep_genes if 'DeprecatedCodeword' not in x]

    # Subset the DataFrame to keep only relevant genes and clusters
    df_sig = df_sig.loc[keep_genes, clusters]

    # drop columns with Nan values
    df_sig = df_sig.dropna(axis=1, how='all')

    df_sig = df_sig.loc[sorted(df_sig.index), sorted(df_sig.columns)]

    # Save the gene expression signatures
    if any(isinstance(dtype, pd.SparseDtype) for dtype in df_sig.dtypes):
        df_sig.sparse.to_dense().to_parquet(os.path.join(
                path_landscape_files,
                f"df_sig{'_' + segmentation_approach if segmentation_approach != 'default' else ''}.parquet"
            ))
    else:
        df_sig.to_parquet(os.path.join(
                path_landscape_files,
                f"df_sig{'_' + segmentation_approach if segmentation_approach != 'default' else ''}.parquet"
            ))

    print("Cluster-specific gene expression signatures saved successfully.")

    return df_sig

create_cluster_and_meta_cluster(technology, path_landscape_files, data_dir=None, segmentation_approach='default')

Creates cell clusters and meta cluster files for visualization. Currently supports only Xenium.

Parameters:

Name	Type	Description	Default
technology	str	The technology used (e.g., "Xenium" or "MERSCOPE"). Currently, only "Xenium" is supported.	required
data_dir	str	Path to the directory containing the Xenium data.	None
path_landscape_files	str	Path to the directory where the cluster and meta cluster files will be saved.	required

Raises:

Type	Description
ValueError	If the specified technology is not supported.
FileNotFoundError	If the required input files are not found.

Source code in src/celldega/pre/__init__.py

def create_cluster_and_meta_cluster(technology, path_landscape_files, data_dir=None, segmentation_approach='default'):
    """
    Creates cell clusters and meta cluster files for visualization.
    Currently supports only Xenium.

    Args:
        technology (str): The technology used (e.g., "Xenium" or "MERSCOPE"). Currently, only "Xenium" is supported.
        data_dir (str): Path to the directory containing the Xenium data.
        path_landscape_files (str): Path to the directory where the cluster and meta cluster files will be saved.

    Raises:
        ValueError: If the specified technology is not supported.
        FileNotFoundError: If the required input files are not found.
    """

    print("\n========Create clusters and meta clusters files========")

    if technology != "Xenium" and technology != 'custom':
        raise ValueError(
            f"Unsupported technology: {technology}. Currently, only 'Xenium' and 'Custom' is supported."
        )

    # Check if the cell metadata file exists
    cell_metadata_path = f"{path_landscape_files}/cell_metadata{'_' + segmentation_approach if segmentation_approach != 'default' else ''}.parquet"

    if not os.path.exists(cell_metadata_path):
        raise FileNotFoundError(
            f"The file '{os.path.basename(cell_metadata_path)}' does not exist in directory '{path_landscape_files}'."
        )

    # Create the cell_clusters directory if it doesn't exist
    cell_clusters_dir = f"{path_landscape_files}/cell_clusters{'_' + segmentation_approach if segmentation_approach != 'default' else ''}"

    if not os.path.exists(cell_clusters_dir):
        os.mkdir(cell_clusters_dir)

    # Load the cell metadata
    meta_cell = pd.read_parquet(cell_metadata_path)

    if technology == 'Xenium':

        # Load the default clustering data (replace this with actual data loading logic)
        default_clustering = pd.read_csv(
            f'{data_dir}/analysis/clustering/gene_expression_graphclust/clusters.csv',
            index_col=0,
        )
        default_clustering.columns = default_clustering.columns.str.lower()

        # Prepare the clustering data
        default_clustering_ini = default_clustering.copy()
        default_clustering_ini['cluster'] = default_clustering_ini['cluster'].astype(
            'string'
        )

        # Align the clustering data with the cell metadata
        default_clustering = pd.DataFrame(index=meta_cell['name'].tolist())
        default_clustering.loc[default_clustering_ini.index.tolist(), 'cluster'] = (
            default_clustering_ini['cluster']
        )

        # Save the clustering data
        default_clustering.to_parquet(f'{cell_clusters_dir}/cluster.parquet')

        # Count the number of cells in each cluster
        ser_counts = default_clustering['cluster'].value_counts()
        clusters = ser_counts.index.tolist()

        # Assign colors to clusters
        palettes = [plt.get_cmap(name).colors for name in plt.colormaps() if "tab" in name]
        flat_colors = [color for palette in palettes for color in palette]
        flat_colors_hex = [to_hex(color) for color in flat_colors]

        colors = [
            (
                flat_colors_hex[i % len(flat_colors_hex)]
                if "Blank" not in cluster
                else "#FFFFFF"
            )
            for i, cluster in enumerate(clusters)
        ]

        # Create the meta cluster DataFrame
        ser_color = pd.Series(colors, index=clusters, name='color')
        meta_cluster = pd.DataFrame(ser_color)
        meta_cluster['count'] = ser_counts

        # Save the meta cluster data
        meta_cluster.to_parquet(f'{cell_clusters_dir}/meta_cluster.parquet')

    if technology == 'custom':

        df_cluster = pd.DataFrame(index=meta_cell['name'].tolist())
        df_cluster['cluster'] = 0
        df_cluster['cluster'] = df_cluster['cluster'].astype('string')
        df_cluster.to_parquet(os.path.join(cell_clusters_dir, "cluster.parquet"))

        meta_cluster = pd.DataFrame(index=['0'])
        meta_cluster.loc['0', 'color'] = '#1f77b4'
        meta_cluster.loc['0', 'count'] = len(meta_cell['name'].tolist())
        meta_cluster.to_parquet(os.path.join(cell_clusters_dir, "meta_cluster.parquet"))

        ser_counts = df_cluster['cluster'].value_counts()
        clusters = ser_counts.index.tolist()

    print("Cell clusters and meta cluster files created successfully.")

    return clusters

create_image_tiles(technology, data_dir, path_landscape_files, image_tile_layer='dapi')

Creates image tiles for visualization from the Xenium morphology image.

Parameters:

Name	Type	Description	Default
technology	str	The technology used (e.g., "Xenium" or "MERSCOPE"). Currently, only "Xenium" is supported.	required
data_dir	str	Path to the directory containing the Xenium data (e.g., morphology_focus_0000.ome.tif).	required
path_landscape_files	str	Path to the directory where the image tiles and pyramid will be saved.	required
image_tile_layer	str	Specifies which image layers to process. Options are 'dapi' (default) or 'all'.	'dapi'

Raises:

Type	Description
ValueError	If the specified technology is not supported or if the image_tile_layer is invalid.
FileNotFoundError	If the required input image file is not found.

Source code in src/celldega/pre/__init__.py

def create_image_tiles(
    technology, data_dir, path_landscape_files, image_tile_layer='dapi'
):
    """
    Creates image tiles for visualization from the Xenium morphology image.

    Args:
        technology (str): The technology used (e.g., "Xenium" or "MERSCOPE"). Currently, only "Xenium" is supported.
        data_dir (str): Path to the directory containing the Xenium data (e.g., morphology_focus_0000.ome.tif).
        path_landscape_files (str): Path to the directory where the image tiles and pyramid will be saved.
        image_tile_layer (str, optional): Specifies which image layers to process. Options are 'dapi' (default) or 'all'.

    Raises:
        ValueError: If the specified technology is not supported or if the image_tile_layer is invalid.
        FileNotFoundError: If the required input image file is not found.
    """

    print("\n========Generating image tiles========")
    if technology != "Xenium":
        raise ValueError(
            f"Unsupported technology: {technology}. Currently, only 'Xenium' is supported."
        )

    if image_tile_layer not in ['dapi', 'all']:
        raise ValueError(
            f"Invalid image_tile_layer: {image_tile_layer}. Must be 'dapi' or 'all'."
        )

    # Define the path to the morphology image
    file_path = f"{data_dir}/morphology_focus/morphology_focus_0000.ome.tif"

    # Check if the morphology image exists
    if not os.path.exists(file_path):
        raise FileNotFoundError(
            f"The file 'morphology_focus_0000.ome.tif' does not exist in directory '{data_dir}'."
        )

    # Process the DAPI channel
    if image_tile_layer == 'dapi' or image_tile_layer == 'all':
        f'generating DAPI image tiles ...'

        if os.path.exists(f'{path_landscape_files}/pyramid_images/dapi_files'):
            pass
        else:

            # Load the morphology image
            img = imread(file_path)

            # Save the DAPI channel to a regular TIFF file
            imsave(f'{path_landscape_files}/dapi_output_regular.tif', img[..., 0])

            # Convert the image to PNG format
            image_png = _convert_to_png(f'{path_landscape_files}/dapi_output_regular.tif')

            # Create a DeepZoom pyramid for the DAPI channel
            make_deepzoom_pyramid(
                image_png,
                f'{path_landscape_files}/pyramid_images/',
                'dapi',
                suffix=".webp[Q=100]",
            )


    # Process additional channels if image_tile_layer is 'all'
    if image_tile_layer == 'all':
        for idx, channel in enumerate(['bound', 'rna', 'prot']):
            print(f'generating {channel} image tiles ...')

            if os.path.exists(f'{path_landscape_files}/pyramid_images/{channel}_files'):
                pass
            else:
                # Extract and process each channel
                image_data = (
                    img[..., idx + 1] * 2
                )  # Adjust intensity for better visualization
                imsave(
                    f'{path_landscape_files}/{channel}_output_regular.tif', image_data
                )

                # Convert the image to PNG format
                image_png = _convert_to_png(f'{path_landscape_files}/{channel}_output_regular.tif')

                # Create a DeepZoom pyramid for the channel
                make_deepzoom_pyramid(
                    image_png,
                    f'{path_landscape_files}/pyramid_images/',
                    channel,
                    suffix=".webp[Q=100]",
                )
    # Remove intermediate files
    intermediate_image_files = glob.glob(f"{path_landscape_files}/*output_regular*")
    if len(intermediate_image_files) != 0: [os.remove(file) for file in intermediate_image_files]

    print("Image tiles created successfully.")

get_image_info(technology, image_tile_layer='dapi')

Retrieve image information for a given technology and image tile layer.

Parameters:

Name	Type	Description	Default
technology	str	The technology for which image information is requested. Currently supports 'Xenium' and 'MERSCOPE'.	required
image_tile_layer	str	The type of image tile layer to retrieve information for. Options are 'dapi' or 'all'. Defaults to 'dapi'.	'dapi'

Returns:

Type	Description
list[dict]	A list of dictionaries containing image information, including name,
list[dict]	button name, and color.

Raises:

Type	Description
ValueError	If the technology is not supported or the image_tile_layer is invalid.

Source code in src/celldega/pre/image_info.py

def get_image_info(technology: str, image_tile_layer: str = 'dapi') -> list[dict]:
    """Retrieve image information for a given technology and image tile layer.

    Args:
        technology: The technology for which image information is requested.
                   Currently supports 'Xenium' and 'MERSCOPE'.
        image_tile_layer: The type of image tile layer to retrieve information for.
                         Options are 'dapi' or 'all'. Defaults to 'dapi'.

    Returns:
        A list of dictionaries containing image information, including name,
        button name, and color.

    Raises:
        ValueError: If the technology is not supported or the image_tile_layer
                   is invalid.
    """
    # Validate technology
    supported_technologies = ['Xenium', 'MERSCOPE']
    if technology not in supported_technologies:
        raise ValueError(
            f"Unsupported technology: {technology}. Supported technologies are: {supported_technologies}."
        )

    # Validate image_tile_layer
    if image_tile_layer not in ['dapi', 'all']:
        raise ValueError(
            f"Invalid image_tile_layer: {image_tile_layer}. Must be 'dapi' or 'all'."
        )

    # Handle 'dapi' case for both Xenium and MERSCOPE
    if image_tile_layer == 'dapi':
        image_info = [{"name": "dapi", "button_name": "DAPI", "color": [0, 0, 255]}]
    # Handle 'all' case (only for Xenium)
    elif image_tile_layer == 'all':
        if technology != 'Xenium':
            raise ValueError(
                f"image_tile_layer='all' is only supported for 'Xenium'. "
                f"Received technology: {technology}."
            )
        image_info = [
            {"name": "dapi", "button_name": "DAPI", "color": [0, 0, 255]},
            {"name": "bound", "button_name": "BOUND", "color": [0, 255, 0]},
            {"name": "rna", "button_name": "RNA", "color": [255, 0, 0]},
            {"name": "prot", "button_name": "PROT", "color": [255, 255, 255]},
        ]

    return image_info

get_max_zoom_level(path_image_pyramid)

Returns the maximum zoom level based on the highest-numbered directory in the specified path.

Parameters:

Name	Type	Description	Default
path_image_pyramid	str	Path to the directory containing zoom level directories.	required

Returns:

Name	Type	Description
int		The maximum zoom level.

Source code in src/celldega/pre/__init__.py

def get_max_zoom_level(path_image_pyramid):
    """Returns the maximum zoom level based on the highest-numbered directory in the specified path.

    Args:
        path_image_pyramid (str): Path to the directory containing zoom level directories.

    Returns:
        int: The maximum zoom level.
    """
    zoom_levels = [
        entry
        for entry in os.listdir(path_image_pyramid)
        if os.path.isdir(os.path.join(path_image_pyramid, entry)) and entry.isdigit()
    ]
    max_pyramid_zoom = max(map(int, zoom_levels)) if zoom_levels else None
    return max_pyramid_zoom

hc(df, filter_N_top=None, norm_col='total', norm_row='zscore')

This function performs hierarchical clustering on the rows and columns of a DataFrame and returns the visualization JSON for the Celldega-Matrix method. This function is developed using the approaches and code adaptations from the Clustergrammer2 project.

Parameters:

Name	Type	Description	Default
df	DataFrame	A Pandas DataFrame.	required
filter_N_top	int	The number of top dimensions to keep after filtering.	None
norm_col	str	The type of normalization to apply to the columns.	'total'
norm_row	str	The type of normalization to apply to the rows.	'zscore'

Returns:

Name	Type	Description
dict		The visualization JSON for the Celldega-Matrix method

Example: See Landscape-Matrix_Xenium notebook.

Source code in src/celldega/clust/__init__.py

def hc(df, filter_N_top=None, norm_col='total', norm_row='zscore'):

  """
  This function performs hierarchical clustering on the rows and columns of a
  DataFrame and returns the visualization JSON for the Celldega-Matrix method.
  This function is developed using the approaches and code adaptations from the
  Clustergrammer2 project.

  Args:
    df (pd.DataFrame): A Pandas DataFrame.
    filter_N_top (int): The number of top dimensions to keep after filtering.
    norm_col (str): The type of normalization to apply to the columns.
    norm_row (str): The type of normalization to apply to the rows.

  Returns:
    dict: The visualization JSON for the Celldega-Matrix method


  Example:
  See [Landscape-Matrix_Xenium](../../../examples/brief_notebooks/Landscape-Matrix_Xenium) notebook.
  """

  net = Network()

  net.load_df(df)

  if filter_N_top is not None:
    net.filter_N_top(axis='row', N_top=5000)

  if norm_col == 'total':
    net.normalize(axis='col', norm_type='umi')

  if norm_row == 'zscore':
    net.normalize(axis='row', norm_type='zscore')

  net.cluster()

  network = net.viz

  return network

main(sample, data_root_dir, tile_size, image_tile_layer, path_landscape_files, use_int_index=True)

Main function to preprocess Xenium data and generate landscape files.

Parameters:

Name	Type	Description	Default
sample	str	Name of the sample (e.g., 'Xenium_V1_human_Pancreas_FFPE_outs').	required
data_root_dir	str	Root directory containing the sample data.	required
tile_size	int	Size of the tiles for transcript and boundary tiles.	required
image_tile_layer	str	Image layers to be tiled. 'dapi' or 'all'.	required
path_landscape_files	str	Directory to save the landscape files.	required

Example

change directory to celldega, and run:

python run_pre_processing.py --sample Xenium_V1_human_Pancreas_FFPE_outs --data_root_dir data --tile_size 250 --image_tile_layer 'dapi' --path_landscape_files notebooks/Xenium_V1_human_Pancreas_FFPE_outs

Source code in src/celldega/pre/run_pre_processing.py

def main(sample, data_root_dir, tile_size, image_tile_layer, path_landscape_files, use_int_index=True):
    """
    Main function to preprocess Xenium data and generate landscape files.

    Args:
        sample (str): Name of the sample (e.g., 'Xenium_V1_human_Pancreas_FFPE_outs').
        data_root_dir (str): Root directory containing the sample data.
        tile_size (int): Size of the tiles for transcript and boundary tiles.
        image_tile_layer (str): Image layers to be tiled. 'dapi' or 'all'.
        path_landscape_files (str): Directory to save the landscape files.

    Example:
        change directory to celldega, and run:

        python run_pre_processing.py \
            --sample Xenium_V1_human_Pancreas_FFPE_outs \
            --data_root_dir data \
            --tile_size 250 \
            --image_tile_layer 'dapi' \
            --path_landscape_files notebooks/Xenium_V1_human_Pancreas_FFPE_outs

    """
    print(f"Starting preprocessing for sample: {sample}")

    # Construct data directory
    data_dir = os.path.join(data_root_dir, sample)

    # Create necessary directories if they don't exist
    for folder in [data_dir, path_landscape_files]:
        if not os.path.exists(folder):
            os.makedirs(folder)
            print(f"Created directory: {folder}")

    # Determine technology based on the presence of experiment.xenium file
    if os.path.exists(os.path.join(data_dir, 'experiment.xenium')):
        technology = 'Xenium'
    else:
        raise ValueError("Unsupported technology. Only Xenium is supported in this script.")

    # Unzip compressed files in Xenium data folder
    dega.pre._xenium_unzipper(data_dir)

    # Check required files for preprocessing
    dega.pre._check_required_files(technology, data_dir)

    # Write transform file
    transformation_matrix = dega.pre.write_xenium_transform(data_dir, path_landscape_files)

    # Make cell image coordinates
    path_transformation_matrix = os.path.join(path_landscape_files, 'micron_to_image_transform.csv')
    path_meta_cell_micron = os.path.join(data_dir, 'cells.csv.gz')
    path_meta_cell_image = os.path.join(path_landscape_files, 'cell_metadata.parquet')
    dega.pre.make_meta_cell_image_coord(
        technology,
        path_transformation_matrix,
        path_meta_cell_micron,
        path_meta_cell_image,
        image_scale=1
    )

    # Calculate CBG
    cbg = dega.pre.read_cbg_mtx(os.path.join(data_dir, 'cell_feature_matrix'))

    # Create cluster-based gene expression
    df_sig = dega.pre.cluster_gene_expression(technology, path_landscape_files, cbg, data_dir)

    # Make meta gene files
    path_output = os.path.join(path_landscape_files, 'meta_gene.parquet')
    dega.pre.make_meta_gene(cbg, path_output)

    # Save CBG gene parquet files
    dega.pre.save_cbg_gene_parquets(path_landscape_files, cbg, verbose=True)

    # Create cluster and meta cluster files
    clusters = dega.pre.create_cluster_and_meta_cluster(technology, path_landscape_files, data_dir)

    # Generate image tiles
    dega.pre.create_image_tiles(technology, data_dir, path_landscape_files, image_tile_layer=image_tile_layer)

    # Generate transcript tiles
    print("\n========Generating transcript tiles========")
    path_trx = os.path.join(data_dir, 'transcripts.parquet')
    path_trx_tiles = os.path.join(path_landscape_files, 'transcript_tiles')
    tile_bounds = dega.pre.make_trx_tiles(
        technology,
        path_trx,
        path_transformation_matrix,
        path_trx_tiles,
        coarse_tile_factor=10,
        tile_size=tile_size,
        chunk_size=100000,
        verbose=False,
        image_scale=1,
        max_workers=2
    )
    print (f"tile bounds: {tile_bounds}")

    # Generate boundary tiles
    print("\n========Generating boundary tiles========")
    path_cell_boundaries = os.path.join(data_dir, 'cell_boundaries.parquet')
    path_output = os.path.join(path_landscape_files, 'cell_segmentation')

    dega.pre.make_cell_boundary_tiles(
            technology,
            path_cell_boundaries,
            path_output,
            path_meta_cell_micron,
            path_transformation_matrix,
            coarse_tile_factor=10,
            tile_size=tile_size,
            tile_bounds=tile_bounds,
            image_scale=1,
            max_workers=2
        )

    # Save landscape parameters
    dega.pre.save_landscape_parameters(
        technology,
        path_landscape_files,
        'dapi_files',
        tile_size=tile_size,
        image_info=dega.pre.get_image_info(technology, image_tile_layer),
        image_format='.webp',
        use_int_index=use_int_index,
    )

    print("Preprocessing completed successfully.")

make_cell_boundary_tiles(technology, path_cell_boundaries, path_output, path_meta_cell_micron=None, path_transformation_matrix=None, coarse_tile_factor=20, tile_size=250, tile_bounds=None, image_scale=1, max_workers=8)

Processes cell boundary data and divides it into spatial tiles based on the provided technology. Reads cell boundary data, applies affine transformations, and divides the data into coarse and fine tiles. The resulting tiles are saved as Parquet files, each containing the geometries of cells in that tile.

Parameters

technology : str The technology used to generate the cell boundary data, e.g., "MERSCOPE", "Xenium", or "custom". path_cell_boundaries : str Path to the file containing the cell boundaries (Parquet format). path_meta_cell_micron : str Path to the file containing cell metadata (CSV format). path_transformation_matrix : str Path to the file containing the transformation matrix (CSV format). path_output : str Directory path where the output files (Parquet files) for each tile will be saved. coarse_tile_factor : int, optional, default=20. scaling factor of each coarse-grain tile comparing to the fine tile size. tile_size : int, optional, default=500 Size of each fine-grain tile in microns. tile_bounds : dict, optional Dictionary containing the minimum and maximum bounds for x and y coordinates. image_scale : float, optional, default=1 Scale factor to apply to the geometry data. max_workers : int, optional, default=8 Maximum number of parallel workers for processing tiles.

Returns

None

Source code in src/celldega/pre/boundary_tile.py

def make_cell_boundary_tiles(
    technology,
    path_cell_boundaries,
    path_output,
    path_meta_cell_micron=None,
    path_transformation_matrix=None,
    coarse_tile_factor=20,
    tile_size=250,
    tile_bounds=None,
    image_scale=1,
    max_workers=8,
):
    """
    Processes cell boundary data and divides it into spatial tiles based on the provided technology.
    Reads cell boundary data, applies affine transformations, and divides the data into coarse and fine tiles.
    The resulting tiles are saved as Parquet files, each containing the geometries of cells in that tile.

    Parameters
    ----------
    technology : str
        The technology used to generate the cell boundary data, e.g., "MERSCOPE", "Xenium", or "custom".
    path_cell_boundaries : str
        Path to the file containing the cell boundaries (Parquet format).
    path_meta_cell_micron : str
        Path to the file containing cell metadata (CSV format).
    path_transformation_matrix : str
        Path to the file containing the transformation matrix (CSV format).
    path_output : str
        Directory path where the output files (Parquet files) for each tile will be saved.
    coarse_tile_factor  : int, optional, default=20.
        scaling factor of each coarse-grain tile comparing to the fine tile size.
    tile_size : int, optional, default=500
        Size of each fine-grain tile in microns.
    tile_bounds : dict, optional
        Dictionary containing the minimum and maximum bounds for x and y coordinates.
    image_scale : float, optional, default=1
        Scale factor to apply to the geometry data.
    max_workers : int, optional, default=8
        Maximum number of parallel workers for processing tiles.

    Returns
    -------
    None
    """

    print("\n========Create cell boundary spatial tiles========")

    # Ensure the output directory exists
    if not os.path.exists(path_output):
        os.makedirs(path_output)

    if technology == 'custom':
        gdf_cells = gpd.read_parquet(path_cell_boundaries)

        # Convert string index to integer index
        cell_str_to_int_mapping = _get_name_mapping(
            path_output.split('/cell_segmentation')[0],  # get the path of landscape files
            layer='boundary',
            segmentation=path_output.split('cell_segmentation_')[1], # get the cell segmentation method, such as cellpose2.
            )
        gdf_cells.index = gdf_cells.index.map(cell_str_to_int_mapping)

        gdf_cells.rename(columns={'geometry_image_space': 'GEOMETRY'}, inplace=True)

        gdf_cells["center_x"] = gdf_cells["GEOMETRY"].apply(lambda geom: geom.centroid.x)
        gdf_cells["center_y"] = gdf_cells["GEOMETRY"].apply(lambda geom: geom.centroid.y)

        transformed_geometries = []

        for geom in gdf_cells['GEOMETRY']:

            if isinstance(geom, Polygon):

                exterior_coords = np.array(geom.exterior.coords)
                transformed_geometries.append([exterior_coords.tolist()])

        gdf_cells["GEOMETRY"] = transformed_geometries

    else:

        transformation_matrix = pd.read_csv(
        path_transformation_matrix, header=None, sep=" "
        ).values

        gdf_cells = get_cell_polygons(
            technology,
            path_cell_boundaries,
            transformation_matrix,
            path_output,
            image_scale,
            path_meta_cell_micron,
        )

        # Convert string index to integer index
        cell_str_to_int_mapping = _get_name_mapping(
            path_transformation_matrix.replace('/micron_to_image_transform.csv',''),
            layer='boundary',
            )
        gdf_cells.index = gdf_cells.index.map(cell_str_to_int_mapping)

        gdf_cells["center_x"] = gdf_cells.geometry.centroid.x
        gdf_cells["center_y"] = gdf_cells.geometry.centroid.y

    # Calculate tile bounds and fine/coarse tiles
    x_min, x_max = tile_bounds["x_min"], tile_bounds["x_max"]
    y_min, y_max = tile_bounds["y_min"], tile_bounds["y_max"]
    n_fine_tiles_x = int(np.ceil((x_max - x_min) / tile_size))
    n_fine_tiles_y = int(np.ceil((y_max - y_min) / tile_size))
    n_coarse_tiles_x = int(np.ceil((x_max - x_min) / (coarse_tile_factor * tile_size)))
    n_coarse_tiles_y = int(np.ceil((y_max - y_min) / (coarse_tile_factor * tile_size)))

    # Process coarse tiles in parallel
    for i in tqdm(range(n_coarse_tiles_x), desc="Processing coarse tiles"):
        coarse_tile_x_min = x_min + i * (coarse_tile_factor * tile_size)
        coarse_tile_x_max = coarse_tile_x_min + (coarse_tile_factor * tile_size)

        for j in range(n_coarse_tiles_y):
            coarse_tile_y_min = y_min + j * (coarse_tile_factor * tile_size)
            coarse_tile_y_max = coarse_tile_y_min + (coarse_tile_factor * tile_size)

            coarse_tile = gdf_cells[
                (gdf_cells["center_x"] >= coarse_tile_x_min)
                & (gdf_cells["center_x"] < coarse_tile_x_max)
                & (gdf_cells["center_y"] >= coarse_tile_y_min)
                & (gdf_cells["center_y"] < coarse_tile_y_max)
            ]
            if not coarse_tile.empty:
                process_fine_boundaries(
                    coarse_tile,
                    i,
                    j,
                    coarse_tile_x_min,
                    coarse_tile_x_max,
                    coarse_tile_y_min,
                    coarse_tile_y_max,
                    tile_size,
                    path_output,
                    x_min,
                    y_min,
                    n_fine_tiles_x,
                    n_fine_tiles_y,
                    max_workers,
                )

    print("Done.")

make_deepzoom_pyramid(image_path, output_path, pyramid_name, tile_size=512, overlap=0, suffix='.jpeg')

Creates a DeepZoom image pyramid from a JPEG image.

Parameters:

Name	Type	Description	Default
image_path	str	Path to the JPEG image file.	required
output_path	str	Directory to save the DeepZoom pyramid.	required
pyramid_name	str	Name of the pyramid directory.	required
tile_size	int	Tile size for the DeepZoom pyramid. Defaults to 512.	512
overlap	int	Overlap size for the DeepZoom pyramid. Defaults to 0.	0
suffix	str	Suffix for the DeepZoom pyramid tiles. Defaults to ".jpeg".	'.jpeg'

Returns:

Type	Description
	None

Source code in src/celldega/pre/__init__.py

def make_deepzoom_pyramid(
    image_path, output_path, pyramid_name, tile_size=512, overlap=0, suffix=".jpeg"
):
    """Creates a DeepZoom image pyramid from a JPEG image.

    Args:
        image_path (str): Path to the JPEG image file.
        output_path (str): Directory to save the DeepZoom pyramid.
        pyramid_name (str): Name of the pyramid directory.
        tile_size (int, optional): Tile size for the DeepZoom pyramid. Defaults to 512.
        overlap (int, optional): Overlap size for the DeepZoom pyramid. Defaults to 0.
        suffix (str, optional): Suffix for the DeepZoom pyramid tiles. Defaults to ".jpeg".

    Returns:
        None
    """
    output_path = Path(output_path)
    image = pyvips.Image.new_from_file(image_path, access="sequential")
    output_path.mkdir(parents=True, exist_ok=True)
    output_path = output_path / pyramid_name
    image.dzsave(output_path, tile_size=tile_size, overlap=overlap, suffix=suffix)

make_meta_cell_image_coord(technology, path_transformation_matrix, path_meta_cell_micron, path_meta_cell_image, image_scale=1)

Applies an affine transformation to cell coordinates in microns and saves the transformed coordinates in pixels.

Parameters

technology : str The technology used to generate the data, Xenium and MERSCOPE are supported. path_transformation_matrix : str Path to the transformation matrix file path_meta_cell_micron : str Path to the meta cell file with coordinates in microns path_meta_cell_image : str Path to save the meta cell file with coordinates in pixels

Returns

None

Examples

make_meta_cell_image_coord( ... technology='Xenium', ... path_transformation_matrix='data/transformation_matrix.csv', ... path_meta_cell_micron='data/meta_cell_micron.csv', ... path_meta_cell_image='data/meta_cell_image.parquet' ... ) Args: technology (str): The technology used to generate the data (e.g., "Xenium" or "MERSCOPE"). path_transformation_matrix (str): Path to the transformation matrix file. path_meta_cell_micron (str): Path to the meta cell file with coordinates in microns. path_meta_cell_image (str): Path to save the meta cell file with coordinates in pixels. image_scale (float): Scaling factor to convert micron coordinates to pixel coordinates.

Returns:

Type	Description
	None

Source code in src/celldega/pre/__init__.py

def make_meta_cell_image_coord(
    technology,
    path_transformation_matrix,
    path_meta_cell_micron,
    path_meta_cell_image,
    image_scale=1
):
    """Applies an affine transformation to cell coordinates in microns and saves the transformed coordinates in pixels.

    Parameters
    ----------
    technology : str
        The technology used to generate the data, Xenium and MERSCOPE are supported.
    path_transformation_matrix : str
        Path to the transformation matrix file
    path_meta_cell_micron : str
        Path to the meta cell file with coordinates in microns
    path_meta_cell_image : str
        Path to save the meta cell file with coordinates in pixels

    Returns
    -------
    None

    Examples
    --------
    >>> make_meta_cell_image_coord(
    ...     technology='Xenium',
    ...     path_transformation_matrix='data/transformation_matrix.csv',
    ...     path_meta_cell_micron='data/meta_cell_micron.csv',
    ...     path_meta_cell_image='data/meta_cell_image.parquet'
    ... )
    Args:
        technology (str): The technology used to generate the data (e.g., "Xenium" or "MERSCOPE").
        path_transformation_matrix (str): Path to the transformation matrix file.
        path_meta_cell_micron (str): Path to the meta cell file with coordinates in microns.
        path_meta_cell_image (str): Path to save the meta cell file with coordinates in pixels.
        image_scale (float): Scaling factor to convert micron coordinates to pixel coordinates.

    Returns:
        None
    """

    print("\n========Make meta cells in pixel space========")
    transformation_matrix = pd.read_csv(
        path_transformation_matrix, header=None, sep=" "
    ).values

    sparse_matrix = csr_matrix(transformation_matrix)

    if technology == "MERSCOPE":
        meta_cell = pd.read_csv(
            path_meta_cell_micron, usecols=["EntityID", "center_x", "center_y"]
        )
        meta_cell = _convert_long_id_to_short(meta_cell)
        meta_cell["name"] = meta_cell["cell_id"]
        meta_cell = meta_cell.set_index('cell_id')

    elif technology == "Xenium":
        usecols = ["cell_id", "x_centroid", "y_centroid"]
        meta_cell = pd.read_csv(path_meta_cell_micron, index_col=0, usecols=usecols)
        meta_cell.columns = ["center_x", "center_y"]
        meta_cell["name"] = pd.Series(meta_cell.index, index=meta_cell.index)

    elif technology == "custom":
        meta_cell = gpd.read_parquet(path_meta_cell_micron)
        meta_cell['center_x'] = meta_cell.centroid.x
        meta_cell['center_y'] = meta_cell.centroid.y
        meta_cell["name"] = pd.Series(meta_cell.index, index=meta_cell.index).astype('str')
        meta_cell.drop(['area', 'centroid'], axis=1, inplace=True)

    # Adding a ones column to accommodate for affine transformation
    meta_cell["ones"] = 1
    points = meta_cell[["center_x", "center_y", "ones"]].values

    # Applying the transformation matrix
    transformed_points = sparse_matrix.dot(points.T).T[:, :2]
    #transformed_points = np.dot(transformation_matrix, points.T).T

    meta_cell["center_x"] = transformed_points[:, 0]
    meta_cell["center_y"] = transformed_points[:, 1]
    meta_cell.drop(columns=["ones"], inplace=True)

    meta_cell["center_x"] = meta_cell["center_x"] / image_scale
    meta_cell["center_y"] = meta_cell["center_y"] / image_scale

    meta_cell["geometry"] = meta_cell.apply(
        lambda row: [row["center_x"], row["center_y"]], axis=1
    )

    if technology == "MERSCOPE":
        meta_cell = meta_cell[["name", "geometry", "EntityID"]]
    else:
        meta_cell = meta_cell[["name", "geometry"]]

    # Check if the 'name' column is unique
    if not meta_cell['name'].is_unique:
        warnings.warn("Duplicate cell names found in meta_cell!", UserWarning)

     # Apply rounding to the GEOMETRY column
    meta_cell['geometry'] = meta_cell['geometry'].apply(_round_nested_coord_list)       

    # Force alphabetically sort by 'name'
    meta_cell = meta_cell.sort_values(by=['name']).reset_index(drop=True)
    meta_cell.to_parquet(path_meta_cell_image, index=False)
    print('Done.')

make_meta_gene(cbg, path_output)

Creates a DataFrame with genes and their assigned colors.

Parameters:

Name	Type	Description	Default
cbg	DataFrame	A sparse DataFrame with genes as columns and barcodes as rows..	required
path_output	str	Path to save the meta gene file.	required

Returns:

Type	Description
	None

Source code in src/celldega/pre/__init__.py

def make_meta_gene(cbg, path_output):
    """Creates a DataFrame with genes and their assigned colors.

    Args:
        cbg (pandas.DataFrame): A sparse DataFrame with genes as columns and barcodes as rows..
        path_output (str): Path to save the meta gene file.

    Returns:
        None
    """

    print("\n========Write meta gene files========")
    genes = cbg.columns.tolist()

    palettes = [plt.get_cmap(name).colors for name in plt.colormaps() if "tab" in name]
    flat_colors = [color for palette in palettes for color in palette]
    flat_colors_hex = [to_hex(color) for color in flat_colors]

    colors = [
        flat_colors_hex[i % len(flat_colors_hex)] if "Blank" not in gene else "#FFFFFF"
        for i, gene in enumerate(genes)
    ]

    ser_color = pd.Series(colors, index=genes)
    meta_gene = calc_meta_gene_data(cbg)
    meta_gene['color'] = ser_color

    sparse_cols = [
        col for col in meta_gene.columns if pd.api.types.is_sparse(meta_gene[col])
    ]
    for col in sparse_cols:
        meta_gene[col] = meta_gene[col].sparse.to_dense()

    # Force alphabetically sort by index
    meta_gene.sort_index(inplace=True)
    meta_gene.to_parquet(path_output)
    print("All meta gene files are succesfully saved.")

make_trx_tiles(technology, path_trx, path_transformation_matrix=None, path_trx_tiles=None, coarse_tile_factor=10, tile_size=250, chunk_size=1000000, verbose=False, image_scale=1, max_workers=1)

Processes transcript data by dividing it into coarse-grain and fine-grain tiles, applying transformations, and saving the results in a parallelized manner.

Parameters

technology : str The technology used for generating the transcript data (e.g., "MERSCOPE" or "Xenium"). path_trx : str Path to the file containing the transcript data. path_transformation_matrix : str Path to the file containing the transformation matrix (CSV file). path_trx_tiles : str Directory path where the output files (Parquet files) for each tile will be saved. coarse_tile_factor : int, optional Scaling factor of each coarse-grain tile comparing to the fine tile size. tile_size : int, optional Size of each fine-grain tile in microns (default is 250). chunk_size : int, optional Number of rows to process per chunk for memory efficiency (default is 1000000). verbose : bool, optional Flag to enable verbose output (default is False). image_scale : float, optional Scale factor to apply to the transcript coordinates (default is 0.5). max_workers : int, optional Maximum number of parallel workers for processing tiles (default is 1).

Returns

dict A dictionary containing the bounds of the processed data in both x and y directions.

Source code in src/celldega/pre/trx_tile.py

def make_trx_tiles(
    technology,
    path_trx,
    path_transformation_matrix=None,
    path_trx_tiles=None,
    coarse_tile_factor=10,
    tile_size=250,
    chunk_size=1000000,
    verbose=False,
    image_scale=1,
    max_workers=1,
):
    """
    Processes transcript data by dividing it into coarse-grain and fine-grain tiles,
    applying transformations, and saving the results in a parallelized manner.

    Parameters
    ----------
    technology : str
        The technology used for generating the transcript data (e.g., "MERSCOPE" or "Xenium").
    path_trx : str
        Path to the file containing the transcript data.
    path_transformation_matrix : str
        Path to the file containing the transformation matrix (CSV file).
    path_trx_tiles : str
        Directory path where the output files (Parquet files) for each tile will be saved.
    coarse_tile_factor : int, optional
        Scaling factor of each coarse-grain tile comparing to the fine tile size.
    tile_size : int, optional
        Size of each fine-grain tile in microns (default is 250).
    chunk_size : int, optional
        Number of rows to process per chunk for memory efficiency (default is 1000000).
    verbose : bool, optional
        Flag to enable verbose output (default is False).
    image_scale : float, optional
        Scale factor to apply to the transcript coordinates (default is 0.5).
    max_workers : int, optional
        Maximum number of parallel workers for processing tiles (default is 1).

    Returns
    -------
    dict
        A dictionary containing the bounds of the processed data in both x and y directions.
    """

    if technology == 'custom':

        x_min, y_min = 0, 0
        x_max, y_max = pl.read_parquet(path_trx).select(
            [
                pl.col("x_image_coords").max().alias("x_max"),
                pl.col("y_image_coords").max().alias("y_max"),
            ]
        ).row(0)

    else:

        if not os.path.exists(path_trx_tiles):
            os.makedirs(path_trx_tiles)

        transformation_matrix = np.loadtxt(path_transformation_matrix)

        # transformed_points = np.dot(points, transformation_matrix.T)[:, :2]
        # sparse_matrix = csr_matrix(transformation_matrix)
        # transformed_points = sparse_matrix.dot(points.T).T[:, :2]

        gene_str_to_int_mapping = _get_name_mapping(
            path_transformation_matrix.replace('/micron_to_image_transform.csv',''),
            layer='transcript',
            )

        trx = transform_transcript_coordinates(
            technology, path_trx, chunk_size, transformation_matrix, image_scale, gene_str_to_int_mapping=gene_str_to_int_mapping,
        )

        # Get min and max x, y values
        x_min, y_min = 0, 0
        x_max, y_max = trx.select(
            [
                pl.col("transformed_x").max().alias("x_max"),
                pl.col("transformed_y").max().alias("y_max"),
            ]
        ).row(0)

        # Calculate the number of fine-grain tiles globally
        n_fine_tiles_x = int(np.ceil((x_max - x_min) / tile_size))
        n_fine_tiles_y = int(np.ceil((y_max - y_min) / tile_size))

        # Calculate the number of coarse-grain tiles
        n_coarse_tiles_x = int(np.ceil((x_max - x_min) / (coarse_tile_factor * tile_size)))
        n_coarse_tiles_y = int(np.ceil((y_max - y_min) / (coarse_tile_factor * tile_size)))

        # Use ThreadPoolExecutor for parallel processing of coarse-grain tiles
        with concurrent.futures.ThreadPoolExecutor(max_workers=max_workers) as executor:
            futures = []
            for i in range(n_coarse_tiles_x):
                coarse_tile_x_min = x_min + i * (coarse_tile_factor * tile_size)
                coarse_tile_x_max = coarse_tile_x_min + (coarse_tile_factor * tile_size)

                for j in range(n_coarse_tiles_y):
                    coarse_tile_y_min = y_min + j * (coarse_tile_factor * tile_size)
                    coarse_tile_y_max = coarse_tile_y_min + (coarse_tile_factor * tile_size)

                    # Submit each coarse tile for parallel processing
                    futures.append(
                        executor.submit(
                            process_coarse_tile,
                            trx,
                            i,
                            j,
                            coarse_tile_x_min,
                            coarse_tile_x_max,
                            coarse_tile_y_min,
                            coarse_tile_y_max,
                            tile_size,
                            path_trx_tiles,
                            x_min,
                            y_min,
                            n_fine_tiles_x,
                            n_fine_tiles_y,
                            max_workers,
                        )
                    )

            # Wait for all coarse tiles to complete
            for future in tqdm(
                concurrent.futures.as_completed(futures),
                desc="Processing coarse tiles",
                unit="tile",
            ):
                future.result()  # Raise exceptions if any occurred during execution

    # Return the tile bounds
    tile_bounds = {
        "x_min": x_min,
        "x_max": x_max,
        "y_min": y_min,
        "y_max": y_max,
    }

    return tile_bounds

numpy_affine_transform(coords, matrix)

Apply affine transformation to numpy coordinates.

Source code in src/celldega/pre/boundary_tile.py

def numpy_affine_transform(coords, matrix):
    """Apply affine transformation to numpy coordinates."""
    # Homogeneous coordinates for affine transformation
    coords = np.hstack([coords, np.ones((coords.shape[0], 1))])
    transformed_coords = coords @ matrix.T
    return transformed_coords[:, :2]  # Drop the homogeneous coordinate

read_cbg_mtx(base_path)

Read the cell-by-gene matrix from the mtx files.

Parameters

base_path : str The base path to the directory containing the mtx files.

Returns

cbg : pandas.DataFrame A sparse DataFrame with genes as columns and barcodes as rows.

Source code in src/celldega/pre/landscape.py

def read_cbg_mtx(base_path):
    """
    Read the cell-by-gene matrix from the mtx files.

    Parameters
    ----------
    base_path : str
        The base path to the directory containing the mtx files.

    Returns
    -------
    cbg : pandas.DataFrame
        A sparse DataFrame with genes as columns and barcodes as rows.
    """
    # print("Reading mtx file from ", base_path)

    # File paths
    barcodes_path = os.path.join(base_path, "barcodes.tsv.gz")
    features_path = os.path.join(base_path, "features.tsv.gz")
    matrix_path = os.path.join(base_path, "matrix.mtx.gz")

    # Read barcodes and features
    barcodes = pd.read_csv(barcodes_path, header=None, compression="gzip")
    features = pd.read_csv(features_path, header=None, compression="gzip", sep="\t")

    # Read the gene expression matrix and transpose it
    # Transpose and convert to CSC format for fast column slicing
    matrix = mmread(matrix_path).transpose().tocsc()

    # Create a sparse DataFrame with genes as columns and barcodes as rows
    cbg = pd.DataFrame.sparse.from_spmatrix(
        matrix, index=barcodes[0], columns=features[1]
    )
    cbg = cbg.rename_axis('__index_level_0__', axis='columns')

    return cbg

save_cbg_gene_parquets(base_path, cbg, verbose=False, segmentation_approach='default')

Save the cell-by-gene matrix as gene-specific Parquet files.

Parameters

base_path : str The base path to the parent directory containing the landscape_files directory. cbg : pandas.DataFrame A sparse DataFrame with genes as columns and barcodes as rows. verbose : bool, optional Whether to print progress information, by default False.

Returns

None

Source code in src/celldega/pre/landscape.py

def save_cbg_gene_parquets(base_path, cbg, verbose=False, segmentation_approach='default'):
    """
    Save the cell-by-gene matrix as gene-specific Parquet files.

    Parameters
    ----------
    base_path : str
        The base path to the parent directory containing the landscape_files directory.
    cbg : pandas.DataFrame
        A sparse DataFrame with genes as columns and barcodes as rows.
    verbose : bool, optional
        Whether to print progress information, by default False.

    Returns
    -------
    None
    """
    output_dir = os.path.join(base_path, f"cbg{f'_{segmentation_approach}' if segmentation_approach !='default' else ''}")
    print(output_dir)
    os.makedirs(output_dir, exist_ok=True)

    # convert cell index from string to integer
    cell_str_to_int_mapping = _get_name_mapping(base_path, layer='boundary', segmentation=segmentation_approach)
    cbg.index = cbg.index.map(cell_str_to_int_mapping)

    for index, gene in enumerate(cbg.columns):
        if verbose and index % 100 == 0:
            print(f"Processing gene {index}: {gene}")

        # Extract the column as a DataFrame as a copy
        col_df = cbg[[gene]].copy()

        # Create a DataFrame necessary to prevent error in to_parquet
        inst_df = pd.DataFrame(
            col_df.values, columns=[gene], index=col_df.index.tolist()
        )

        # Replace 0 with NA and drop rows where all values are NA
        inst_df.replace(0, pd.NA, inplace=True)
        inst_df.dropna(how="all", inplace=True)

        # Save to Parquet if DataFrame is not empty
        if not inst_df.empty:
            output_path = os.path.join(output_dir, f"{gene}.parquet")
            inst_df.to_parquet(output_path)

    print("All gene-specific parquet files are succesfully saved.")

save_landscape_parameters(technology, path_landscape_files, image_name='dapi_files', tile_size=1000, image_info={}, image_format='.webp', use_int_index=False, segmentation_approach='default')

Saves the landscape parameters to a JSON file.

Parameters:

Name	Type	Description	Default
technology	str	The technology used to generate the data.	required
path_landscape_files	str	Path to the directory where landscape files are stored.	required
image_name	str	Name of the image directory. Defaults to "dapi_files".	'dapi_files'
tile_size	int	Tile size for the image pyramid. Defaults to 1000.	1000
image_info	dict	Additional image metadata. Defaults to {}.	{}
image_format	str	Format of the image files. Defaults to ".webp".	'.webp'
use_int_index	bool	Use integer name for cell_tile and trx_tile.	False

Returns:

Type	Description
	None

Source code in src/celldega/pre/__init__.py

def save_landscape_parameters(
    technology,
    path_landscape_files,
    image_name="dapi_files",
    tile_size=1000,
    image_info={},
    image_format='.webp',
    use_int_index=False,
    segmentation_approach='default'
):
    """Saves the landscape parameters to a JSON file.

    Args:
        technology (str): The technology used to generate the data.
        path_landscape_files (str): Path to the directory where landscape files are stored.
        image_name (str, optional): Name of the image directory. Defaults to "dapi_files".
        tile_size (int, optional): Tile size for the image pyramid. Defaults to 1000.
        image_info (dict, optional): Additional image metadata. Defaults to {}.
        image_format (str, optional): Format of the image files. Defaults to ".webp".
        use_int_index (bool, optional): Use integer name for cell_tile and trx_tile.

    Returns:
        None
    """

    print("\n========Save landscape parameters========")
    path_image_pyramid = f"{path_landscape_files}/pyramid_images/{image_name}"
    max_pyramid_zoom = get_max_zoom_level(path_image_pyramid)

    path_landscape_parameters = f"{path_landscape_files}/landscape_parameters.json"

    if technology != 'custom':

        landscape_parameters = {
                "technology": technology,
                "segmentation_approach": [segmentation_approach],
                "max_pyramid_zoom": max_pyramid_zoom,
                "tile_size": tile_size,
                "image_info": image_info,
                "image_format": image_format,
                "use_int_index":use_int_index,
            }

    else:

        with open(path_landscape_parameters, "r") as file:
            landscape_parameters = json.load(file)

        landscape_parameters['segmentation_approach'].append(segmentation_approach)

    with open(path_landscape_parameters, "w") as file:
        json.dump(landscape_parameters, file, indent=4)

    print('Done.')

write_xenium_transform(data_dir, path_landscape_files, transform_fname='micron_to_image_transform.csv')

Extracts the transformation matrix from the Xenium cells.zarr.zip file and saves it as a CSV file.

Parameters:

Name	Type	Description	Default
data_dir	str	Path to the directory containing the Xenium data (e.g., cells.zarr.zip).	required
path_landscape_files	str	Path to the directory where the transformation matrix CSV will be saved.	required
transform_fname	str	Name of the output CSV file. Defaults to "micron_to_image_transform.csv".	'micron_to_image_transform.csv'

Returns:

Type	Description
	numpy.ndarray: The full transformation matrix extracted from the Xenium cells.zarr.zip file.

Raises:

Type	Description
FileNotFoundError	If the cells.zarr.zip file does not exist in the specified data_dir.
KeyError	If the transformation matrix is not found in the Zarr file under the expected path.
Exception	If an unexpected error occurs while processing the Zarr file.

Source code in src/celldega/pre/__init__.py

def write_xenium_transform(
    data_dir, path_landscape_files, transform_fname="micron_to_image_transform.csv"
):
    """
    Extracts the transformation matrix from the Xenium cells.zarr.zip file and saves it as a CSV file.

    Args:
        data_dir (str): Path to the directory containing the Xenium data (e.g., cells.zarr.zip).
        path_landscape_files (str): Path to the directory where the transformation matrix CSV will be saved.
        transform_fname (str, optional): Name of the output CSV file. Defaults to "micron_to_image_transform.csv".

    Returns:
        numpy.ndarray: The full transformation matrix extracted from the Xenium cells.zarr.zip file.

    Raises:
        FileNotFoundError: If the cells.zarr.zip file does not exist in the specified `data_dir`.
        KeyError: If the transformation matrix is not found in the Zarr file under the expected path.
        Exception: If an unexpected error occurs while processing the Zarr file.
    """

    print("\n========Write xenium transform file from the Zarr folder========")
    # Path to the cells.zarr.zip file
    cells_zarr_path = os.path.join(data_dir, "cells.zarr.zip")

    # Check if the cells.zarr.zip file exists
    if not os.path.exists(cells_zarr_path):
        raise FileNotFoundError(
            f"The file 'cells.zarr.zip' does not exist in directory '{data_dir}'."
        )

    # Function to open a Zarr file
    def open_zarr(path: str) -> zarr.Group:
        store = (
            zarr.ZipStore(path, mode="r")
            if path.endswith(".zip")
            else zarr.DirectoryStore(path)
        )
        return zarr.group(store=store)

    try:
        # Open the cells Zarr file
        root = open_zarr(cells_zarr_path)

        # Extract the transformation matrix
        transformation_matrix = root['masks']['homogeneous_transform'][:]

        # Save the transformation matrix as a CSV file
        output_path = os.path.join(path_landscape_files, transform_fname)
        pd.DataFrame(transformation_matrix[:3, :3]).to_csv(
            output_path, sep=" ", header=False, index=False
        )

        print(f"Transformation matrix saved to '{output_path}'.")
    except KeyError as e:
        raise KeyError(
            f"Could not find the transformation matrix in the Zarr file: {e}"
        )
    except Exception as e:
        raise Exception(f"An error occurred while processing the Zarr file: {e}")

    return transformation_matrix