Direct links to results

TF-MoDISco results

Summary of motifs

TOMTOM matches to motifs

Sample of seqlets for each motif

import os
import sys
sys.path.append(os.path.abspath("/users/amtseng/tfmodisco/src/"))
from tfmodisco.run_tfmodisco import import_shap_scores, import_tfmodisco_results
from motif.read_motifs import pfm_info_content, pfm_to_pwm, trim_motif_by_ic
from motif.match_motifs import match_motifs_to_database
from util import figure_to_vdom_image, import_peak_table, import_profiles
import plot.viz_sequence as viz_sequence
import numpy as np
import h5py
import pandas as pd
import sklearn.cluster
import scipy.cluster.hierarchy
import matplotlib.pyplot as plt
import matplotlib.font_manager as font_manager
import vdom.helpers as vdomh
from IPython.display import display
import tqdm
tqdm.tqdm_notebook()

/users/amtseng/miniconda3/envs/tfmodisco-mini/lib/python3.7/site-packages/ipykernel_launcher.py:19: TqdmDeprecationWarning: This function will be removed in tqdm==5.0.0
Please use `tqdm.notebook.tqdm` instead of `tqdm.tqdm_notebook`

<tqdm.notebook.tqdm_notebook at 0x7f5ca36be0d0>

# Plotting defaults
font_manager.fontManager.ttflist.extend(
    font_manager.createFontList(
        font_manager.findSystemFonts(fontpaths="/users/amtseng/modules/fonts")
    )
)
plot_params = {
    "figure.titlesize": 22,
    "axes.titlesize": 22,
    "axes.labelsize": 20,
    "legend.fontsize": 18,
    "xtick.labelsize": 16,
    "ytick.labelsize": 16,
    "font.family": "Roboto",
    "font.weight": "bold"
}
plt.rcParams.update(plot_params)

/users/amtseng/miniconda3/envs/tfmodisco-mini/lib/python3.7/site-packages/ipykernel_launcher.py:4: MatplotlibDeprecationWarning: 
The createFontList function was deprecated in Matplotlib 3.2 and will be removed two minor releases later. Use FontManager.addfont instead.
  after removing the cwd from sys.path.

Define constants and paths

# Define parameters/fetch arguments
tf_name = os.environ["TFM_TF_NAME"]
preds_path = os.environ["TFM_PRED_PATH"]
shap_scores_path = os.environ["TFM_SHAP_PATH"]
tfm_results_path = os.environ["TFM_TFM_PATH"]
if "TFM_TASK_INDEX" in os.environ:
    task_index = int(os.environ["TFM_TASK_INDEX"])
else:
    task_index = None
hyp_score_key = os.environ["TFM_HYP_SCORE_KEY"]
if "TFM_MOTIF_CACHE" in os.environ:
    tfm_motifs_cache_dir = os.environ["TFM_MOTIF_CACHE"]
else:
    tfm_motifs_cache_dir = None

print("TF name: %s" % tf_name)
print("Predictions/performance path: %s" % preds_path)
print("DeepSHAP scores path: %s" % shap_scores_path)
print("TF-MoDISco results path: %s" % tfm_results_path)
print("Task index: %s" % task_index)
print("Importance score key: %s" % hyp_score_key)
print("Saved TF-MoDISco-derived motifs cache: %s" % tfm_motifs_cache_dir)

TF name: MAX
Predictions/performance path: /users/amtseng/tfmodisco/results/peak_predictions/multitask_profile_finetune/MAX_multitask_profile_finetune_fold1/MAX_multitask_profile_finetune_fold1_pred_perf.h5
DeepSHAP scores path: /users/amtseng/tfmodisco/results/importance_scores/multitask_profile_finetune/MAX_multitask_profile_finetune_fold1/MAX_multitask_profile_finetune_fold1_imp_scores.h5
TF-MoDISco results path: /users/amtseng/tfmodisco/results/tfmodisco/multitask_profile_finetune/MAX_multitask_profile_finetune_fold1/MAX_multitask_profile_finetune_fold1_count_tfm.h5
Task index: None
Importance score key: count_hyp_scores
Saved TF-MoDISco-derived motifs cache: /users/amtseng/tfmodisco/results/reports/tfmodisco_results//cache/multitask_profile_finetune/MAX_multitask_profile_finetune_fold1/MAX_multitask_profile_finetune_fold1_count

# Define constants
input_length, profile_length = 2114, 1000
shap_score_center_size = 400
profile_display_center_size = 400

base_path = "/users/amtseng/tfmodisco/"
data_path = os.path.join(base_path, "data/processed/ENCODE/")
labels_path = os.path.join(data_path, "labels/%s" % tf_name)

# Paths to original called peaks
all_peak_beds = sorted([item for item in os.listdir(labels_path) if item.endswith(".bed.gz")])
if task_index is None:
    peak_bed_paths = [os.path.join(labels_path, item) for item in all_peak_beds]
else:
    peak_bed_paths = [os.path.join(labels_path, all_peak_beds[task_index])]

if tfm_motifs_cache_dir:
    os.makedirs(tfm_motifs_cache_dir, exist_ok=True)

Helper functions

def extract_profiles_and_coords(
    seqlets_arr, one_hot_seqs, hyp_scores, true_profs, pred_profs, pred_coords,
    input_length, profile_length, input_center_cut_size, profile_center_cut_size,
    task_index=None
):
    """
    From the seqlets object of a TF-MoDISco pattern's seqlets and alignments,
    extracts the predicted and observed profiles of the model, as well as the
    set of coordinates for the seqlets.
    Arguments:
        `seqlets_arr`: a TF-MoDISco pattern's seqlets object array (N-array)
        `one_hot_seqs`: an N x R x 4 array of input sequences, where R is
            the cut centered size
        `hyp_scores`: an N x R x 4 array of hypothetical importance scores
        `true_profs`: an N x T x O x 2 array of true profile counts
        `pred_profs`: an N x T x O x 2 array of predicted profile probabilities
        `pred_coords`: an N x 3 object array of coordinates for the input sequence
            underlying the predictions
        `input_length`: length of original input sequences, I
        `profile_length`: length of profile predictions, O
        `input_center_cut_size`: centered cut size of SHAP scores used
        `profile_center_cut_size`: size to cut profiles to when returning them, P
        `task_index`: index of task to focus on for profiles; if None, returns
            profiles for all tasks
    Returns an N x (T or 1) x P x 2 array of true profile counts, an
    N x (T or 1) x P x 2 array of predicted profile probabilities, an N x Q x 4
    array of one-hot seqlet sequences, an N x Q x 4 array of hypothetical seqlet
    importance scores, and an N x 3 object array of seqlet coordinates, where P
    is the profile cut size and Q is the seqlet length. Returned profiles are
    centered at the same center as the seqlets.
    Note that it is important that the seqlet indices match exactly with the indices
    out of the N. This should be the exact sequences in the original SHAP scores.
    """
    true_seqlet_profs, pred_seqlet_profs, seqlet_seqs, seqlet_hyps, seqlet_coords = [], [], [], [], []
    
    def seqlet_coord_to_profile_coord(seqlet_coord):
        return seqlet_coord + ((input_length - input_center_cut_size) // 2) - ((input_length - profile_length) // 2)
    
    def seqlet_coord_to_input_coord(seqlet_coord):
        return seqlet_coord + ((input_length - input_center_cut_size) // 2)
        
    # For each seqlet, fetch the true/predicted profiles
    for seqlet in seqlets_arr:
        coord_index = seqlet.coor.example_idx
        seqlet_start = seqlet.coor.start
        seqlet_end = seqlet.coor.end
        seqlet_rc = seqlet.coor.is_revcomp
        
        # Get indices of profile to cut out
        seqlet_center = (seqlet_start + seqlet_end) // 2
        prof_center = seqlet_coord_to_profile_coord(seqlet_center)
        prof_start = prof_center - (profile_center_cut_size // 2)
        prof_end = prof_start + profile_center_cut_size
        
        if task_index is None or true_profs.shape[1] == 1:
            # Use all tasks if the predictions only have 1 task to begin with
            task_start, task_end = None, None
        else:
            task_start, task_end = task_index, task_index + 1
            
        true_prof = true_profs[coord_index, task_start:task_end, prof_start:prof_end]  # (T or 1) x P x 2
        pred_prof = pred_profs[coord_index, task_start:task_end, prof_start:prof_end]  # (T or 1) x P x 2
        
        true_seqlet_profs.append(true_prof)
        pred_seqlet_profs.append(pred_prof)
        
        # The one-hot-sequences and hypothetical scores are assumed to already by cut/centered,
        # so the indices match the seqlet indices
        if seqlet_rc:
            seqlet_seqs.append(np.flip(one_hot_seqs[coord_index, seqlet_start:seqlet_end], axis=(0, 1)))
            seqlet_hyps.append(np.flip(hyp_scores[coord_index, seqlet_start:seqlet_end], axis=(0, 1)))
        else:
            seqlet_seqs.append(one_hot_seqs[coord_index, seqlet_start:seqlet_end])
            seqlet_hyps.append(hyp_scores[coord_index, seqlet_start:seqlet_end])
            
        # Get the coordinates of the seqlet based on the input coordinates
        inp_start = seqlet_coord_to_input_coord(seqlet_start)
        inp_end = seqlet_coord_to_input_coord(seqlet_end)
        chrom, start, _ = pred_coords[coord_index]
        seqlet_coords.append([chrom, start + inp_start, start + inp_end])
    
    return np.stack(true_seqlet_profs), np.stack(pred_seqlet_profs), np.stack(seqlet_seqs), np.stack(seqlet_hyps), np.array(seqlet_coords, dtype=object)

def plot_profiles(seqlet_true_profs, seqlet_pred_profs, kmeans_clusters=5, save_path=None):
    """
    Plots the given profiles with a heatmap.
    Arguments:
        `seqlet_true_profs`: an N x O x 2 NumPy array of true profiles, either as raw
            counts or probabilities (they will be normalized)
        `seqlet_pred_profs`: an N x O x 2 NumPy array of predicted profiles, either as
            raw counts or probabilities (they will be normalized)
        `kmeans_cluster`: when displaying profile heatmaps, there will be this
            many clusters
        `save_path`: if provided, save the profile matrices here
    Returns the figure.
    """
    assert len(seqlet_true_profs.shape) == 3
    assert seqlet_true_profs.shape == seqlet_pred_profs.shape
    num_profs, width, _ = seqlet_true_profs.shape

    # First, normalize the profiles along the output profile dimension
    def normalize(arr, axis=0):
        arr_sum = np.sum(arr, axis=axis, keepdims=True)
        arr_sum[arr_sum == 0] = 1  # If 0, keep 0 as the quotient instead of dividing by 0
        return arr / arr_sum
    true_profs_norm = normalize(seqlet_true_profs, axis=1)
    pred_profs_norm = normalize(seqlet_pred_profs, axis=1)

    # Compute the mean profiles across all examples
    true_profs_mean = np.mean(true_profs_norm, axis=0)
    pred_profs_mean = np.mean(pred_profs_norm, axis=0)

    # Perform k-means clustering on the predicted profiles, with the strands pooled
    kmeans_clusters = max(5, num_profs // 50)  # Set number of clusters based on number of profiles, with minimum
    kmeans = sklearn.cluster.KMeans(n_clusters=kmeans_clusters)
    cluster_assignments = kmeans.fit_predict(
        np.reshape(pred_profs_norm, (pred_profs_norm.shape[0], -1))
    )

    # Perform hierarchical clustering on the cluster centers to determine optimal ordering
    kmeans_centers = kmeans.cluster_centers_
    cluster_order = scipy.cluster.hierarchy.leaves_list(
        scipy.cluster.hierarchy.optimal_leaf_ordering(
            scipy.cluster.hierarchy.linkage(kmeans_centers, method="centroid"), kmeans_centers
        )
    )

    # Order the profiles so that the cluster assignments follow the optimal ordering
    cluster_inds = []
    for cluster_id in cluster_order:
        cluster_inds.append(np.where(cluster_assignments == cluster_id)[0])
    cluster_inds = np.concatenate(cluster_inds)

    # Compute a matrix of profiles, normalized to the maximum height, ordered by clusters
    def make_profile_matrix(flat_profs, order_inds):
        matrix = flat_profs[order_inds]
        maxes = np.max(matrix, axis=1, keepdims=True)
        maxes[maxes == 0] = 1  # If 0, keep 0 as the quotient instead of dividing by 0
        return matrix / maxes
    true_matrix = make_profile_matrix(true_profs_norm, cluster_inds)
    pred_matrix = make_profile_matrix(pred_profs_norm, cluster_inds)
    
    if save_path:
        np.savez_compressed(
            true_profs_mean=true_profs_mean, pred_profs_mean=pred_profs_mean,
            true_matrix=true_matrix, pred_matrix=pred_matrix
        )

    # Create a figure with the right dimensions
    mean_height = 4
    heatmap_height = min(num_profs * 0.004, 8)
    fig_height = mean_height + (2 * heatmap_height)
    fig, ax = plt.subplots(
        3, 2, figsize=(16, fig_height), sharex=True,
        gridspec_kw={
            "width_ratios": [1, 1],
            "height_ratios": [mean_height / fig_height, heatmap_height / fig_height, heatmap_height / fig_height]
        }
    )

    # Plot the average predictions
    ax[0, 0].plot(true_profs_mean[:, 0], color="darkslateblue")
    ax[0, 0].plot(-true_profs_mean[:, 1], color="darkorange")
    ax[0, 1].plot(pred_profs_mean[:, 0], color="darkslateblue")
    ax[0, 1].plot(-pred_profs_mean[:, 1], color="darkorange")

    # Set axes on average predictions
    max_mean_val = max(np.max(true_profs_mean), np.max(pred_profs_mean))
    mean_ylim = max_mean_val * 1.05  # Make 5% higher
    ax[0, 0].set_title("True profiles")
    ax[0, 0].set_ylabel("Average probability")
    ax[0, 1].set_title("Predicted profiles")
    for j in (0, 1):
        ax[0, j].set_ylim(-mean_ylim, mean_ylim)
        ax[0, j].label_outer()

    # Plot the heatmaps
    ax[1, 0].imshow(true_matrix[:, :, 0], interpolation="nearest", aspect="auto", cmap="Blues")
    ax[1, 1].imshow(pred_matrix[:, :, 0], interpolation="nearest", aspect="auto", cmap="Blues")
    ax[2, 0].imshow(true_matrix[:, :, 1], interpolation="nearest", aspect="auto", cmap="Oranges")
    ax[2, 1].imshow(pred_matrix[:, :, 1], interpolation="nearest", aspect="auto", cmap="Oranges")

    # Set axes on heatmaps
    for i in (1, 2):
        for j in (0, 1):
            ax[i, j].set_yticks([])
            ax[i, j].set_yticklabels([])
            ax[i, j].label_outer()
    width = true_matrix.shape[1]
    delta = 100
    num_deltas = (width // 2) // delta
    labels = list(range(max(-width // 2, -num_deltas * delta), min(width // 2, num_deltas * delta) + 1, delta))
    tick_locs = [label + max(width // 2, num_deltas * delta) for label in labels]
    for j in (0, 1):
        ax[2, j].set_xticks(tick_locs)
        ax[2, j].set_xticklabels(labels)
        ax[2, j].set_xlabel("Distance from seqlet center (bp)")

    fig.tight_layout()
    plt.show()
    return fig

def get_summit_distances(coords, peak_table):
    """
    Given a set of coordinates, computes the distance of the center of each
    coordinate to the nearest summit.
    Arguments:
        `coords`: an N x 3 object array of coordinates
        `peak_table`: a 10-column table of peak data, as imported by
            `import_peak_table`
    Returns and N-array of integers, which is the distance of each coordinate
    midpoint to the nearest coordinate.
    """
    chroms = coords[:, 0]
    midpoints = (coords[:, 1] + coords[:, 2]) // 2
    dists = []
    for i in range(len(coords)):
        chrom = chroms[i]
        midpoint = midpoints[i]
        rows = peak_table[peak_table["chrom"] == chrom]
        dist_arr = (midpoint - rows["summit"]).values
        min_dist = dist_arr[np.argmin(np.abs(dist_arr))]
        dists.append(min_dist)
    return np.array(dists)

def plot_summit_dists(summit_dists):
    """
    Plots the distribution of seqlet distances to summits.
    Arguments:
        `summit_dists`: the array of distances as returned by
            `get_summit_distances`
    Returns the figure.
    """
    fig = plt.figure(figsize=(8, 6))
    num_bins = max(len(summit_dists) // 30, 20)
    plt.hist(summit_dists, bins=num_bins, color="purple")
    plt.title("Histogram of distance of seqlets to peak summits")
    plt.xlabel("Signed distance from seqlet center to nearest peak summit (bp)")
    plt.show()
    return fig

Import SHAP scores, profile predictions, and TF-MoDISco results

# Import SHAP coordinates and one-hot sequences
hyp_scores, _, one_hot_seqs, shap_coords = import_shap_scores(shap_scores_path, hyp_score_key, center_cut_size=shap_score_center_size)
# This cuts the sequences/scores off just as how TF-MoDISco saw them, but the coordinates are uncut

Importing SHAP scores:   0%|          | 0/204 [00:00<?, ?it/s]
Importing SHAP scores:   0%|          | 1/204 [00:01<04:08,  1.22s/it]
Importing SHAP scores:   1%|          | 2/204 [00:02<03:51,  1.15s/it]
Importing SHAP scores:   1%|▏         | 3/204 [00:02<02:21,  1.42it/s]
Importing SHAP scores:   2%|▏         | 4/204 [00:03<03:02,  1.09it/s]
Importing SHAP scores:   2%|▏         | 5/204 [00:04<03:12,  1.04it/s]
Importing SHAP scores:   3%|▎         | 6/204 [00:04<02:17,  1.44it/s]
Importing SHAP scores:   3%|▎         | 7/204 [00:06<02:44,  1.20it/s]
Importing SHAP scores:   4%|▍         | 8/204 [00:07<02:57,  1.10it/s]
Importing SHAP scores:   4%|▍         | 9/204 [00:07<02:11,  1.49it/s]
Importing SHAP scores:   5%|▍         | 10/204 [00:08<02:38,  1.22it/s]
Importing SHAP scores:   5%|▌         | 11/204 [00:08<01:59,  1.62it/s]
Importing SHAP scores:   6%|▌         | 12/204 [00:09<02:35,  1.23it/s]
Importing SHAP scores:   6%|▋         | 13/204 [00:10<02:49,  1.13it/s]
Importing SHAP scores:   7%|▋         | 14/204 [00:11<02:07,  1.50it/s]
Importing SHAP scores:   7%|▋         | 15/204 [00:11<02:10,  1.45it/s]
Importing SHAP scores:   8%|▊         | 16/204 [00:12<02:10,  1.44it/s]
Importing SHAP scores:   8%|▊         | 17/204 [00:12<01:39,  1.88it/s]
Importing SHAP scores:   9%|▉         | 18/204 [00:13<01:44,  1.78it/s]
Importing SHAP scores:   9%|▉         | 19/204 [00:13<01:21,  2.27it/s]
Importing SHAP scores:  10%|▉         | 20/204 [00:14<01:33,  1.98it/s]
Importing SHAP scores:  10%|█         | 21/204 [00:14<01:48,  1.69it/s]
Importing SHAP scores:  11%|█         | 22/204 [00:15<01:23,  2.17it/s]
Importing SHAP scores:  11%|█▏        | 23/204 [00:15<01:42,  1.76it/s]
Importing SHAP scores:  12%|█▏        | 24/204 [00:16<01:58,  1.52it/s]
Importing SHAP scores:  12%|█▏        | 25/204 [00:16<01:30,  1.97it/s]
Importing SHAP scores:  13%|█▎        | 26/204 [00:17<01:37,  1.83it/s]
Importing SHAP scores:  13%|█▎        | 27/204 [00:18<01:43,  1.71it/s]
Importing SHAP scores:  14%|█▎        | 28/204 [00:18<01:20,  2.19it/s]
Importing SHAP scores:  14%|█▍        | 29/204 [00:19<01:32,  1.88it/s]
Importing SHAP scores:  15%|█▍        | 30/204 [00:19<01:13,  2.37it/s]
Importing SHAP scores:  15%|█▌        | 31/204 [00:19<01:25,  2.03it/s]
Importing SHAP scores:  16%|█▌        | 32/204 [00:21<01:56,  1.48it/s]
Importing SHAP scores:  16%|█▌        | 33/204 [00:21<01:29,  1.91it/s]
Importing SHAP scores:  17%|█▋        | 34/204 [00:22<01:54,  1.48it/s]
Importing SHAP scores:  17%|█▋        | 35/204 [00:23<02:21,  1.19it/s]
Importing SHAP scores:  18%|█▊        | 36/204 [00:23<01:52,  1.49it/s]
Importing SHAP scores:  18%|█▊        | 37/204 [00:24<02:15,  1.23it/s]
Importing SHAP scores:  19%|█▊        | 38/204 [00:25<01:42,  1.61it/s]
Importing SHAP scores:  19%|█▉        | 39/204 [00:26<02:15,  1.21it/s]
Importing SHAP scores:  20%|█▉        | 40/204 [00:27<02:37,  1.04it/s]
Importing SHAP scores:  20%|██        | 41/204 [00:27<01:57,  1.38it/s]
Importing SHAP scores:  21%|██        | 42/204 [00:28<02:19,  1.16it/s]
Importing SHAP scores:  21%|██        | 43/204 [00:30<02:30,  1.07it/s]
Importing SHAP scores:  22%|██▏       | 44/204 [00:30<01:54,  1.39it/s]
Importing SHAP scores:  22%|██▏       | 45/204 [00:31<02:17,  1.16it/s]
Importing SHAP scores:  23%|██▎       | 46/204 [00:31<01:42,  1.53it/s]
Importing SHAP scores:  23%|██▎       | 47/204 [00:32<02:03,  1.27it/s]
Importing SHAP scores:  24%|██▎       | 48/204 [00:33<02:18,  1.13it/s]
Importing SHAP scores:  24%|██▍       | 49/204 [00:34<01:43,  1.49it/s]
Importing SHAP scores:  25%|██▍       | 50/204 [00:35<02:05,  1.23it/s]
Importing SHAP scores:  25%|██▌       | 51/204 [00:35<01:57,  1.30it/s]
Importing SHAP scores:  25%|██▌       | 52/204 [00:36<01:30,  1.68it/s]
Importing SHAP scores:  26%|██▌       | 53/204 [00:36<01:36,  1.56it/s]
Importing SHAP scores:  26%|██▋       | 54/204 [00:37<01:40,  1.49it/s]
Importing SHAP scores:  27%|██▋       | 55/204 [00:37<01:17,  1.93it/s]
Importing SHAP scores:  27%|██▋       | 56/204 [00:38<01:25,  1.73it/s]
Importing SHAP scores:  28%|██▊       | 57/204 [00:38<01:07,  2.19it/s]
Importing SHAP scores:  28%|██▊       | 58/204 [00:39<01:20,  1.82it/s]
Importing SHAP scores:  29%|██▉       | 59/204 [00:40<01:27,  1.66it/s]
Importing SHAP scores:  29%|██▉       | 60/204 [00:40<01:08,  2.12it/s]
Importing SHAP scores:  30%|██▉       | 61/204 [00:41<01:21,  1.76it/s]
Importing SHAP scores:  30%|███       | 62/204 [00:41<01:26,  1.65it/s]
Importing SHAP scores:  31%|███       | 63/204 [00:41<01:07,  2.09it/s]
Importing SHAP scores:  31%|███▏      | 64/204 [00:42<01:26,  1.61it/s]
Importing SHAP scores:  32%|███▏      | 65/204 [00:43<01:07,  2.06it/s]
Importing SHAP scores:  32%|███▏      | 66/204 [00:43<01:16,  1.80it/s]
Importing SHAP scores:  33%|███▎      | 67/204 [00:44<01:23,  1.63it/s]
Importing SHAP scores:  33%|███▎      | 68/204 [00:44<01:05,  2.09it/s]
Importing SHAP scores:  34%|███▍      | 69/204 [00:45<01:14,  1.82it/s]
Importing SHAP scores:  34%|███▍      | 70/204 [00:46<01:20,  1.66it/s]
Importing SHAP scores:  35%|███▍      | 71/204 [00:46<01:02,  2.12it/s]
Importing SHAP scores:  35%|███▌      | 72/204 [00:47<01:14,  1.77it/s]
Importing SHAP scores:  36%|███▌      | 73/204 [00:47<00:58,  2.23it/s]
Importing SHAP scores:  36%|███▋      | 74/204 [00:48<01:11,  1.83it/s]
Importing SHAP scores:  37%|███▋      | 75/204 [00:48<01:17,  1.66it/s]
Importing SHAP scores:  37%|███▋      | 76/204 [00:49<01:03,  2.00it/s]
Importing SHAP scores:  38%|███▊      | 77/204 [00:49<01:11,  1.77it/s]
Importing SHAP scores:  38%|███▊      | 78/204 [00:50<01:14,  1.69it/s]
Importing SHAP scores:  39%|███▊      | 79/204 [00:50<00:58,  2.15it/s]
Importing SHAP scores:  39%|███▉      | 80/204 [00:51<01:05,  1.88it/s]
Importing SHAP scores:  40%|███▉      | 81/204 [00:52<01:17,  1.58it/s]
Importing SHAP scores:  40%|████      | 82/204 [00:52<01:01,  1.98it/s]
Importing SHAP scores:  41%|████      | 83/204 [00:53<01:08,  1.76it/s]
Importing SHAP scores:  41%|████      | 84/204 [00:53<00:53,  2.23it/s]
Importing SHAP scores:  42%|████▏     | 85/204 [00:53<01:01,  1.92it/s]
Importing SHAP scores:  42%|████▏     | 86/204 [00:54<01:08,  1.72it/s]
Importing SHAP scores:  43%|████▎     | 87/204 [00:54<00:53,  2.20it/s]
Importing SHAP scores:  43%|████▎     | 88/204 [00:55<00:59,  1.95it/s]
Importing SHAP scores:  44%|████▎     | 89/204 [00:56<01:07,  1.70it/s]
Importing SHAP scores:  44%|████▍     | 90/204 [00:56<00:52,  2.18it/s]
Importing SHAP scores:  45%|████▍     | 91/204 [00:57<01:00,  1.87it/s]
Importing SHAP scores:  45%|████▌     | 92/204 [00:57<00:49,  2.28it/s]
Importing SHAP scores:  46%|████▌     | 93/204 [00:58<00:59,  1.88it/s]
Importing SHAP scores:  46%|████▌     | 94/204 [00:58<01:02,  1.75it/s]
Importing SHAP scores:  47%|████▋     | 95/204 [00:58<00:48,  2.23it/s]
Importing SHAP scores:  47%|████▋     | 96/204 [00:59<00:58,  1.86it/s]
Importing SHAP scores:  48%|████▊     | 97/204 [01:00<01:03,  1.68it/s]
Importing SHAP scores:  48%|████▊     | 98/204 [01:00<00:49,  2.13it/s]
Importing SHAP scores:  49%|████▊     | 99/204 [01:01<00:55,  1.88it/s]
Importing SHAP scores:  49%|████▉     | 100/204 [01:01<00:58,  1.77it/s]
Importing SHAP scores:  50%|████▉     | 101/204 [01:01<00:46,  2.23it/s]
Importing SHAP scores:  50%|█████     | 102/204 [01:02<00:59,  1.71it/s]
Importing SHAP scores:  50%|█████     | 103/204 [01:03<00:46,  2.18it/s]
Importing SHAP scores:  51%|█████     | 104/204 [01:03<00:55,  1.81it/s]
Importing SHAP scores:  51%|█████▏    | 105/204 [01:04<00:59,  1.65it/s]
Importing SHAP scores:  52%|█████▏    | 106/204 [01:04<00:47,  2.07it/s]
Importing SHAP scores:  52%|█████▏    | 107/204 [01:05<00:51,  1.90it/s]
Importing SHAP scores:  53%|█████▎    | 108/204 [01:06<00:55,  1.72it/s]
Importing SHAP scores:  53%|█████▎    | 109/204 [01:06<00:44,  2.14it/s]
Importing SHAP scores:  54%|█████▍    | 110/204 [01:06<00:50,  1.85it/s]
Importing SHAP scores:  54%|█████▍    | 111/204 [01:07<00:40,  2.31it/s]
Importing SHAP scores:  55%|█████▍    | 112/204 [01:08<00:52,  1.76it/s]
Importing SHAP scores:  55%|█████▌    | 113/204 [01:08<00:59,  1.53it/s]
Importing SHAP scores:  56%|█████▌    | 114/204 [01:09<00:47,  1.90it/s]
Importing SHAP scores:  56%|█████▋    | 115/204 [01:09<00:52,  1.69it/s]
Importing SHAP scores:  57%|█████▋    | 116/204 [01:10<00:54,  1.62it/s]
Importing SHAP scores:  57%|█████▋    | 117/204 [01:10<00:41,  2.09it/s]
Importing SHAP scores:  58%|█████▊    | 118/204 [01:11<00:47,  1.82it/s]
Importing SHAP scores:  58%|█████▊    | 119/204 [01:11<00:36,  2.30it/s]
Importing SHAP scores:  59%|█████▉    | 120/204 [01:12<00:41,  2.01it/s]
Importing SHAP scores:  59%|█████▉    | 121/204 [01:13<00:48,  1.73it/s]
Importing SHAP scores:  60%|█████▉    | 122/204 [01:13<00:37,  2.20it/s]
Importing SHAP scores:  60%|██████    | 123/204 [01:13<00:44,  1.82it/s]
Importing SHAP scores:  61%|██████    | 124/204 [01:14<00:47,  1.69it/s]
Importing SHAP scores:  61%|██████▏   | 125/204 [01:14<00:36,  2.14it/s]
Importing SHAP scores:  62%|██████▏   | 126/204 [01:15<00:43,  1.80it/s]
Importing SHAP scores:  62%|██████▏   | 127/204 [01:16<00:46,  1.65it/s]
Importing SHAP scores:  63%|██████▎   | 128/204 [01:16<00:36,  2.11it/s]
Importing SHAP scores:  63%|██████▎   | 129/204 [01:17<00:40,  1.86it/s]
Importing SHAP scores:  64%|██████▎   | 130/204 [01:17<00:31,  2.34it/s]
Importing SHAP scores:  64%|██████▍   | 131/204 [01:18<00:37,  1.94it/s]
Importing SHAP scores:  65%|██████▍   | 132/204 [01:18<00:45,  1.60it/s]
Importing SHAP scores:  65%|██████▌   | 133/204 [01:19<00:36,  1.94it/s]
Importing SHAP scores:  66%|██████▌   | 134/204 [01:19<00:40,  1.75it/s]
Importing SHAP scores:  66%|██████▌   | 135/204 [01:20<00:41,  1.65it/s]
Importing SHAP scores:  67%|██████▋   | 136/204 [01:20<00:32,  2.09it/s]
Importing SHAP scores:  67%|██████▋   | 137/204 [01:21<00:44,  1.49it/s]
Importing SHAP scores:  68%|██████▊   | 138/204 [01:22<00:34,  1.93it/s]
Importing SHAP scores:  68%|██████▊   | 139/204 [01:22<00:36,  1.80it/s]
Importing SHAP scores:  69%|██████▊   | 140/204 [01:23<00:37,  1.69it/s]
Importing SHAP scores:  69%|██████▉   | 141/204 [01:23<00:30,  2.06it/s]
Importing SHAP scores:  70%|██████▉   | 142/204 [01:24<00:34,  1.78it/s]
Importing SHAP scores:  70%|███████   | 143/204 [01:25<00:37,  1.64it/s]
Importing SHAP scores:  71%|███████   | 144/204 [01:25<00:28,  2.10it/s]
Importing SHAP scores:  71%|███████   | 145/204 [01:25<00:31,  1.85it/s]
Importing SHAP scores:  72%|███████▏  | 146/204 [01:26<00:24,  2.34it/s]
Importing SHAP scores:  72%|███████▏  | 147/204 [01:26<00:29,  1.92it/s]
Importing SHAP scores:  73%|███████▎  | 148/204 [01:27<00:31,  1.76it/s]
Importing SHAP scores:  73%|███████▎  | 149/204 [01:27<00:24,  2.24it/s]
Importing SHAP scores:  74%|███████▎  | 150/204 [01:28<00:29,  1.85it/s]
Importing SHAP scores:  74%|███████▍  | 151/204 [01:29<00:31,  1.69it/s]
Importing SHAP scores:  75%|███████▍  | 152/204 [01:29<00:24,  2.11it/s]
Importing SHAP scores:  75%|███████▌  | 153/204 [01:30<00:27,  1.83it/s]
Importing SHAP scores:  75%|███████▌  | 154/204 [01:30<00:29,  1.71it/s]
Importing SHAP scores:  76%|███████▌  | 155/204 [01:30<00:23,  2.08it/s]
Importing SHAP scores:  76%|███████▋  | 156/204 [01:31<00:26,  1.81it/s]
Importing SHAP scores:  77%|███████▋  | 157/204 [01:31<00:20,  2.30it/s]
Importing SHAP scores:  77%|███████▋  | 158/204 [01:32<00:23,  1.92it/s]
Importing SHAP scores:  78%|███████▊  | 159/204 [01:33<00:26,  1.67it/s]
Importing SHAP scores:  78%|███████▊  | 160/204 [01:33<00:20,  2.13it/s]
Importing SHAP scores:  79%|███████▉  | 161/204 [01:34<00:26,  1.64it/s]
Importing SHAP scores:  79%|███████▉  | 162/204 [01:35<00:26,  1.57it/s]
Importing SHAP scores:  80%|███████▉  | 163/204 [01:35<00:20,  2.02it/s]
Importing SHAP scores:  80%|████████  | 164/204 [01:35<00:21,  1.84it/s]
Importing SHAP scores:  81%|████████  | 165/204 [01:36<00:16,  2.33it/s]
Importing SHAP scores:  81%|████████▏ | 166/204 [01:36<00:19,  1.96it/s]
Importing SHAP scores:  82%|████████▏ | 167/204 [01:37<00:20,  1.79it/s]
Importing SHAP scores:  82%|████████▏ | 168/204 [01:37<00:15,  2.27it/s]
Importing SHAP scores:  83%|████████▎ | 169/204 [01:38<00:18,  1.92it/s]
Importing SHAP scores:  83%|████████▎ | 170/204 [01:39<00:20,  1.64it/s]
Importing SHAP scores:  84%|████████▍ | 171/204 [01:39<00:15,  2.09it/s]
Importing SHAP scores:  84%|████████▍ | 172/204 [01:40<00:17,  1.81it/s]
Importing SHAP scores:  85%|████████▍ | 173/204 [01:40<00:18,  1.70it/s]
Importing SHAP scores:  85%|████████▌ | 174/204 [01:40<00:13,  2.15it/s]
Importing SHAP scores:  86%|████████▌ | 175/204 [01:41<00:15,  1.87it/s]
Importing SHAP scores:  86%|████████▋ | 176/204 [01:41<00:11,  2.35it/s]
Importing SHAP scores:  87%|████████▋ | 177/204 [01:42<00:13,  2.00it/s]
Importing SHAP scores:  87%|████████▋ | 178/204 [01:43<00:14,  1.86it/s]
Importing SHAP scores:  88%|████████▊ | 179/204 [01:43<00:10,  2.29it/s]
Importing SHAP scores:  88%|████████▊ | 180/204 [01:44<00:13,  1.72it/s]
Importing SHAP scores:  89%|████████▊ | 181/204 [01:44<00:14,  1.60it/s]
Importing SHAP scores:  89%|████████▉ | 182/204 [01:45<00:10,  2.06it/s]
Importing SHAP scores:  90%|████████▉ | 183/204 [01:45<00:11,  1.88it/s]
Importing SHAP scores:  90%|█████████ | 184/204 [01:45<00:08,  2.36it/s]
Importing SHAP scores:  91%|█████████ | 185/204 [01:46<00:09,  1.99it/s]
Importing SHAP scores:  91%|█████████ | 186/204 [01:47<00:11,  1.60it/s]
Importing SHAP scores:  92%|█████████▏| 187/204 [01:47<00:08,  2.06it/s]
Importing SHAP scores:  92%|█████████▏| 188/204 [01:48<00:10,  1.59it/s]
Importing SHAP scores:  93%|█████████▎| 189/204 [01:49<00:10,  1.48it/s]
Importing SHAP scores:  93%|█████████▎| 190/204 [01:49<00:07,  1.89it/s]
Importing SHAP scores:  94%|█████████▎| 191/204 [01:50<00:07,  1.73it/s]
Importing SHAP scores:  94%|█████████▍| 192/204 [01:50<00:05,  2.21it/s]
Importing SHAP scores:  95%|█████████▍| 193/204 [01:51<00:06,  1.58it/s]
Importing SHAP scores:  95%|█████████▌| 194/204 [01:52<00:06,  1.49it/s]
Importing SHAP scores:  96%|█████████▌| 195/204 [01:52<00:04,  1.92it/s]
Importing SHAP scores:  96%|█████████▌| 196/204 [01:53<00:04,  1.77it/s]
Importing SHAP scores:  97%|█████████▋| 197/204 [01:53<00:04,  1.63it/s]
Importing SHAP scores:  97%|█████████▋| 198/204 [01:54<00:02,  2.07it/s]
Importing SHAP scores:  98%|█████████▊| 199/204 [01:54<00:02,  1.79it/s]
Importing SHAP scores:  98%|█████████▊| 200/204 [01:55<00:02,  1.71it/s]
Importing SHAP scores:  99%|█████████▊| 201/204 [01:55<00:01,  2.19it/s]
Importing SHAP scores:  99%|█████████▉| 202/204 [01:56<00:01,  1.91it/s]
Importing SHAP scores: 100%|█████████▉| 203/204 [01:56<00:00,  2.40it/s]
Importing SHAP scores: 100%|██████████| 204/204 [01:56<00:00,  1.75it/s]

# Import the set of all profiles and their coordinates
true_profs, pred_profs, all_pred_coords = import_profiles(preds_path)

Importing predictions:   0%|          | 0/204 [00:00<?, ?it/s]
Importing predictions:   0%|          | 1/204 [00:04<14:24,  4.26s/it]
Importing predictions:   1%|          | 2/204 [00:07<12:28,  3.70s/it]
Importing predictions:   1%|▏         | 3/204 [00:10<11:42,  3.50s/it]
Importing predictions:   2%|▏         | 4/204 [00:17<16:07,  4.84s/it]
Importing predictions:   2%|▏         | 5/204 [00:20<14:06,  4.25s/it]
Importing predictions:   3%|▎         | 6/204 [00:24<12:50,  3.89s/it]
Importing predictions:   3%|▎         | 7/204 [00:30<15:54,  4.84s/it]
Importing predictions:   4%|▍         | 8/204 [00:34<14:06,  4.32s/it]
Importing predictions:   4%|▍         | 9/204 [00:37<12:53,  3.97s/it]
Importing predictions:   5%|▍         | 10/204 [00:44<15:40,  4.85s/it]
Importing predictions:   5%|▌         | 11/204 [00:47<13:58,  4.34s/it]
Importing predictions:   6%|▌         | 12/204 [00:50<12:46,  3.99s/it]
Importing predictions:   6%|▋         | 13/204 [00:57<15:19,  4.81s/it]
Importing predictions:   7%|▋         | 14/204 [01:00<13:41,  4.32s/it]
Importing predictions:   7%|▋         | 15/204 [01:03<12:33,  3.98s/it]
Importing predictions:   8%|▊         | 16/204 [01:10<15:10,  4.84s/it]
Importing predictions:   8%|▊         | 17/204 [01:13<13:32,  4.35s/it]
Importing predictions:   9%|▉         | 18/204 [01:16<12:23,  4.00s/it]
Importing predictions:   9%|▉         | 19/204 [01:20<11:35,  3.76s/it]
Importing predictions:  10%|▉         | 20/204 [01:27<14:36,  4.76s/it]
Importing predictions:  10%|█         | 21/204 [01:30<13:06,  4.30s/it]
Importing predictions:  11%|█         | 22/204 [01:33<12:03,  3.97s/it]
Importing predictions:  11%|█▏        | 23/204 [01:40<14:57,  4.96s/it]
Importing predictions:  12%|█▏        | 24/204 [01:44<13:18,  4.44s/it]
Importing predictions:  12%|█▏        | 25/204 [01:47<12:08,  4.07s/it]
Importing predictions:  13%|█▎        | 26/204 [01:56<16:22,  5.52s/it]
Importing predictions:  13%|█▎        | 27/204 [01:59<14:16,  4.84s/it]
Importing predictions:  14%|█▎        | 28/204 [02:02<12:47,  4.36s/it]
Importing predictions:  14%|█▍        | 29/204 [02:11<16:21,  5.61s/it]
Importing predictions:  15%|█▍        | 30/204 [02:14<14:16,  4.92s/it]
Importing predictions:  15%|█▌        | 31/204 [02:17<12:54,  4.48s/it]
Importing predictions:  16%|█▌        | 32/204 [02:25<15:05,  5.26s/it]
Importing predictions:  16%|█▌        | 33/204 [02:28<13:15,  4.65s/it]
Importing predictions:  17%|█▋        | 34/204 [02:31<11:59,  4.23s/it]
Importing predictions:  17%|█▋        | 35/204 [02:38<14:22,  5.10s/it]
Importing predictions:  18%|█▊        | 36/204 [02:42<12:51,  4.59s/it]
Importing predictions:  18%|█▊        | 37/204 [02:45<11:44,  4.22s/it]
Importing predictions:  19%|█▊        | 38/204 [02:48<10:52,  3.93s/it]
Importing predictions:  19%|█▉        | 39/204 [02:56<14:21,  5.22s/it]
Importing predictions:  20%|█▉        | 40/204 [03:00<12:36,  4.61s/it]
Importing predictions:  20%|██        | 41/204 [03:03<11:24,  4.20s/it]
Importing predictions:  21%|██        | 42/204 [03:10<13:27,  4.98s/it]
Importing predictions:  21%|██        | 43/204 [03:13<11:56,  4.45s/it]
Importing predictions:  22%|██▏       | 44/204 [03:16<10:53,  4.09s/it]
Importing predictions:  22%|██▏       | 45/204 [03:23<12:48,  4.83s/it]
Importing predictions:  23%|██▎       | 46/204 [03:26<11:27,  4.35s/it]
Importing predictions:  23%|██▎       | 47/204 [03:29<10:28,  4.00s/it]
Importing predictions:  24%|██▎       | 48/204 [03:37<13:06,  5.04s/it]
Importing predictions:  24%|██▍       | 49/204 [03:40<11:36,  4.50s/it]
Importing predictions:  25%|██▍       | 50/204 [03:43<10:32,  4.11s/it]
Importing predictions:  25%|██▌       | 51/204 [03:50<12:59,  5.09s/it]
Importing predictions:  25%|██▌       | 52/204 [03:54<11:28,  4.53s/it]
Importing predictions:  26%|██▌       | 53/204 [03:57<10:30,  4.18s/it]
Importing predictions:  26%|██▋       | 54/204 [04:11<18:01,  7.21s/it]
Importing predictions:  27%|██▋       | 55/204 [04:18<17:49,  7.18s/it]
Importing predictions:  27%|██▋       | 56/204 [04:26<18:11,  7.38s/it]
Importing predictions:  28%|██▊       | 57/204 [04:35<19:25,  7.93s/it]
Importing predictions:  28%|██▊       | 58/204 [04:51<24:33, 10.09s/it]
Importing predictions:  29%|██▉       | 59/204 [04:59<23:16,  9.63s/it]
Importing predictions:  29%|██▉       | 60/204 [05:06<21:23,  8.91s/it]
Importing predictions:  30%|██▉       | 61/204 [05:22<25:53, 10.87s/it]
Importing predictions:  30%|███       | 62/204 [05:29<23:18,  9.85s/it]
Importing predictions:  31%|███       | 63/204 [05:39<23:14,  9.89s/it]
Importing predictions:  31%|███▏      | 64/204 [05:53<25:59, 11.14s/it]
Importing predictions:  32%|███▏      | 65/204 [06:03<25:08, 10.85s/it]
Importing predictions:  32%|███▏      | 66/204 [06:11<22:41,  9.86s/it]
Importing predictions:  33%|███▎      | 67/204 [06:28<27:16, 11.95s/it]
Importing predictions:  33%|███▎      | 68/204 [06:38<26:09, 11.54s/it]
Importing predictions:  34%|███▍      | 69/204 [06:47<23:52, 10.61s/it]
Importing predictions:  34%|███▍      | 70/204 [07:06<29:19, 13.13s/it]
Importing predictions:  35%|███▍      | 71/204 [07:10<23:17, 10.51s/it]
Importing predictions:  35%|███▌      | 72/204 [07:16<19:47,  8.99s/it]
Importing predictions:  36%|███▌      | 73/204 [07:26<20:11,  9.25s/it]
Importing predictions:  36%|███▋      | 74/204 [07:38<22:07, 10.22s/it]
Importing predictions:  37%|███▋      | 75/204 [07:43<18:26,  8.57s/it]
Importing predictions:  37%|███▋      | 76/204 [07:49<17:06,  8.02s/it]
Importing predictions:  38%|███▊      | 77/204 [07:58<17:30,  8.27s/it]
Importing predictions:  38%|███▊      | 78/204 [08:03<15:01,  7.16s/it]
Importing predictions:  39%|███▊      | 79/204 [08:07<13:10,  6.33s/it]
Importing predictions:  39%|███▉      | 80/204 [08:17<15:25,  7.46s/it]
Importing predictions:  40%|███▉      | 81/204 [08:22<13:14,  6.46s/it]
Importing predictions:  40%|████      | 82/204 [08:25<11:10,  5.49s/it]
Importing predictions:  41%|████      | 83/204 [08:31<11:48,  5.85s/it]
Importing predictions:  41%|████      | 84/204 [08:35<10:08,  5.07s/it]
Importing predictions:  42%|████▏     | 85/204 [08:38<08:57,  4.52s/it]
Importing predictions:  42%|████▏     | 86/204 [08:45<10:14,  5.21s/it]
Importing predictions:  43%|████▎     | 87/204 [08:48<08:59,  4.61s/it]
Importing predictions:  43%|████▎     | 88/204 [08:51<08:06,  4.20s/it]
Importing predictions:  44%|████▎     | 89/204 [08:59<10:20,  5.40s/it]
Importing predictions:  44%|████▍     | 90/204 [09:03<08:59,  4.74s/it]
Importing predictions:  45%|████▍     | 91/204 [09:06<08:03,  4.27s/it]
Importing predictions:  45%|████▌     | 92/204 [09:09<07:22,  3.95s/it]
Importing predictions:  46%|████▌     | 93/204 [09:17<09:43,  5.26s/it]
Importing predictions:  46%|████▌     | 94/204 [09:21<08:39,  4.72s/it]
Importing predictions:  47%|████▋     | 95/204 [09:24<07:45,  4.27s/it]
Importing predictions:  47%|████▋     | 96/204 [09:32<09:32,  5.30s/it]
Importing predictions:  48%|████▊     | 97/204 [09:35<08:20,  4.68s/it]
Importing predictions:  48%|████▊     | 98/204 [09:38<07:28,  4.23s/it]
Importing predictions:  49%|████▊     | 99/204 [09:46<09:07,  5.21s/it]
Importing predictions:  49%|████▉     | 100/204 [09:49<08:02,  4.64s/it]
Importing predictions:  50%|████▉     | 101/204 [09:52<07:14,  4.22s/it]
Importing predictions:  50%|█████     | 102/204 [10:00<09:17,  5.46s/it]
Importing predictions:  50%|█████     | 103/204 [10:04<08:18,  4.94s/it]
Importing predictions:  51%|█████     | 104/204 [10:08<07:25,  4.46s/it]
Importing predictions:  51%|█████▏    | 105/204 [10:15<08:54,  5.40s/it]
Importing predictions:  52%|█████▏    | 106/204 [10:18<07:44,  4.74s/it]
Importing predictions:  52%|█████▏    | 107/204 [10:22<06:58,  4.32s/it]
Importing predictions:  53%|█████▎    | 108/204 [10:29<08:15,  5.16s/it]
Importing predictions:  53%|█████▎    | 109/204 [10:32<07:14,  4.57s/it]
Importing predictions:  54%|█████▍    | 110/204 [10:35<06:31,  4.16s/it]
Importing predictions:  54%|█████▍    | 111/204 [10:38<06:02,  3.90s/it]
Importing predictions:  55%|█████▍    | 112/204 [10:45<07:24,  4.83s/it]
Importing predictions:  55%|█████▌    | 113/204 [10:49<06:35,  4.34s/it]
Importing predictions:  56%|█████▌    | 114/204 [10:52<06:01,  4.02s/it]
Importing predictions:  56%|█████▋    | 115/204 [10:59<07:19,  4.94s/it]
Importing predictions:  57%|█████▋    | 116/204 [11:02<06:29,  4.42s/it]
Importing predictions:  57%|█████▋    | 117/204 [11:05<05:52,  4.06s/it]
Importing predictions:  58%|█████▊    | 118/204 [11:12<07:04,  4.94s/it]
Importing predictions:  58%|█████▊    | 119/204 [11:16<06:16,  4.43s/it]
Importing predictions:  59%|█████▉    | 120/204 [11:19<05:41,  4.06s/it]
Importing predictions:  59%|█████▉    | 121/204 [11:26<07:00,  5.07s/it]
Importing predictions:  60%|█████▉    | 122/204 [11:30<06:10,  4.52s/it]
Importing predictions:  60%|██████    | 123/204 [11:33<05:34,  4.13s/it]
Importing predictions:  61%|██████    | 124/204 [11:40<06:37,  4.96s/it]
Importing predictions:  61%|██████▏   | 125/204 [11:43<05:51,  4.44s/it]
Importing predictions:  62%|██████▏   | 126/204 [11:46<05:18,  4.09s/it]
Importing predictions:  62%|██████▏   | 127/204 [11:53<06:27,  5.04s/it]
Importing predictions:  63%|██████▎   | 128/204 [11:57<05:41,  4.50s/it]
Importing predictions:  63%|██████▎   | 129/204 [12:00<05:08,  4.11s/it]
Importing predictions:  64%|██████▎   | 130/204 [12:03<04:45,  3.86s/it]
Importing predictions:  64%|██████▍   | 131/204 [12:11<06:02,  4.97s/it]
Importing predictions:  65%|██████▍   | 132/204 [12:14<05:20,  4.45s/it]
Importing predictions:  65%|██████▌   | 133/204 [12:17<04:49,  4.08s/it]
Importing predictions:  66%|██████▌   | 134/204 [12:25<06:05,  5.22s/it]
Importing predictions:  66%|██████▌   | 135/204 [12:28<05:19,  4.62s/it]
Importing predictions:  67%|██████▋   | 136/204 [12:32<04:45,  4.20s/it]
Importing predictions:  67%|██████▋   | 137/204 [12:39<05:47,  5.19s/it]
Importing predictions:  68%|██████▊   | 138/204 [12:42<05:04,  4.61s/it]
Importing predictions:  68%|██████▊   | 139/204 [12:45<04:32,  4.19s/it]
Importing predictions:  69%|██████▊   | 140/204 [12:53<05:23,  5.06s/it]
Importing predictions:  69%|██████▉   | 141/204 [12:56<04:44,  4.51s/it]
Importing predictions:  70%|██████▉   | 142/204 [12:59<04:15,  4.13s/it]
Importing predictions:  70%|███████   | 143/204 [13:06<05:00,  4.92s/it]
Importing predictions:  71%|███████   | 144/204 [13:09<04:24,  4.41s/it]
Importing predictions:  71%|███████   | 145/204 [13:12<03:59,  4.07s/it]
Importing predictions:  72%|███████▏  | 146/204 [13:16<03:41,  3.81s/it]
Importing predictions:  72%|███████▏  | 147/204 [13:24<04:52,  5.14s/it]
Importing predictions:  73%|███████▎  | 148/204 [13:27<04:15,  4.57s/it]
Importing predictions:  73%|███████▎  | 149/204 [13:30<03:49,  4.16s/it]
Importing predictions:  74%|███████▎  | 150/204 [13:37<04:27,  4.96s/it]
Importing predictions:  74%|███████▍  | 151/204 [13:40<03:55,  4.45s/it]
Importing predictions:  75%|███████▍  | 152/204 [13:43<03:32,  4.08s/it]
Importing predictions:  75%|███████▌  | 153/204 [13:51<04:15,  5.02s/it]
Importing predictions:  75%|███████▌  | 154/204 [13:54<03:44,  4.49s/it]
Importing predictions:  76%|███████▌  | 155/204 [13:57<03:21,  4.11s/it]
Importing predictions:  76%|███████▋  | 156/204 [14:04<03:58,  4.97s/it]
Importing predictions:  77%|███████▋  | 157/204 [14:07<03:29,  4.45s/it]
Importing predictions:  77%|███████▋  | 158/204 [14:11<03:08,  4.09s/it]
Importing predictions:  78%|███████▊  | 159/204 [14:18<03:43,  4.97s/it]
Importing predictions:  78%|███████▊  | 160/204 [14:21<03:15,  4.45s/it]
Importing predictions:  79%|███████▉  | 161/204 [14:24<02:55,  4.08s/it]
Importing predictions:  79%|███████▉  | 162/204 [14:31<03:29,  4.98s/it]
Importing predictions:  80%|███████▉  | 163/204 [14:34<03:02,  4.45s/it]
Importing predictions:  80%|████████  | 164/204 [14:38<02:43,  4.09s/it]
Importing predictions:  81%|████████  | 165/204 [14:41<02:29,  3.83s/it]
Importing predictions:  81%|████████▏ | 166/204 [14:49<03:16,  5.18s/it]
Importing predictions:  82%|████████▏ | 167/204 [14:52<02:50,  4.60s/it]
Importing predictions:  82%|████████▏ | 168/204 [14:56<02:30,  4.19s/it]
Importing predictions:  83%|████████▎ | 169/204 [15:03<03:01,  5.18s/it]
Importing predictions:  83%|████████▎ | 170/204 [15:06<02:36,  4.59s/it]
Importing predictions:  84%|████████▍ | 171/204 [15:10<02:18,  4.18s/it]
Importing predictions:  84%|████████▍ | 172/204 [15:18<02:53,  5.42s/it]
Importing predictions:  85%|████████▍ | 173/204 [15:21<02:27,  4.76s/it]
Importing predictions:  85%|████████▌ | 174/204 [15:24<02:08,  4.30s/it]
Importing predictions:  86%|████████▌ | 175/204 [15:32<02:37,  5.43s/it]
Importing predictions:  86%|████████▋ | 176/204 [15:36<02:13,  4.78s/it]
Importing predictions:  87%|████████▋ | 177/204 [15:39<01:56,  4.31s/it]
Importing predictions:  87%|████████▋ | 178/204 [15:47<02:19,  5.35s/it]
Importing predictions:  88%|████████▊ | 179/204 [15:50<01:57,  4.71s/it]
Importing predictions:  88%|████████▊ | 180/204 [15:53<01:42,  4.27s/it]
Importing predictions:  89%|████████▊ | 181/204 [16:02<02:07,  5.56s/it]
Importing predictions:  89%|████████▉ | 182/204 [16:05<01:46,  4.86s/it]
Importing predictions:  90%|████████▉ | 183/204 [16:08<01:31,  4.37s/it]
Importing predictions:  90%|█████████ | 184/204 [16:11<01:20,  4.03s/it]
Importing predictions:  91%|█████████ | 185/204 [16:19<01:37,  5.14s/it]
Importing predictions:  91%|█████████ | 186/204 [16:22<01:22,  4.56s/it]
Importing predictions:  92%|█████████▏| 187/204 [16:26<01:10,  4.17s/it]
Importing predictions:  92%|█████████▏| 188/204 [16:33<01:23,  5.24s/it]
Importing predictions:  93%|█████████▎| 189/204 [16:37<01:09,  4.64s/it]
Importing predictions:  93%|█████████▎| 190/204 [16:40<00:59,  4.22s/it]
Importing predictions:  94%|█████████▎| 191/204 [16:47<01:06,  5.10s/it]
Importing predictions:  94%|█████████▍| 192/204 [16:50<00:54,  4.54s/it]
Importing predictions:  95%|█████████▍| 193/204 [16:53<00:45,  4.15s/it]
Importing predictions:  95%|█████████▌| 194/204 [17:02<00:54,  5.43s/it]
Importing predictions:  96%|█████████▌| 195/204 [17:05<00:42,  4.77s/it]
Importing predictions:  96%|█████████▌| 196/204 [17:08<00:34,  4.31s/it]
Importing predictions:  97%|█████████▋| 197/204 [17:15<00:35,  5.12s/it]
Importing predictions:  97%|█████████▋| 198/204 [17:19<00:27,  4.57s/it]
Importing predictions:  98%|█████████▊| 199/204 [17:22<00:20,  4.16s/it]
Importing predictions:  98%|█████████▊| 200/204 [17:29<00:20,  5.08s/it]
Importing predictions:  99%|█████████▊| 201/204 [17:32<00:13,  4.52s/it]
Importing predictions:  99%|█████████▉| 202/204 [17:36<00:08,  4.13s/it]
Importing predictions: 100%|█████████▉| 203/204 [17:39<00:03,  3.86s/it]
Importing predictions: 100%|██████████| 204/204 [17:42<00:00,  5.21s/it]

# Import the set of peaks
peak_table = import_peak_table(peak_bed_paths)

# Subset the predicted profiles/coordinates to the task-specific SHAP coordinates/scores
shap_coords_table = pd.DataFrame(shap_coords, columns=["chrom", "start", "end"])
pred_coords_table = pd.DataFrame(all_pred_coords, columns=["chrom", "start", "end"])

subset_inds = pred_coords_table.reset_index().drop_duplicates(["chrom", "start", "end"]).merge(
    shap_coords_table.reset_index(), on=["chrom", "start", "end"]
).sort_values("index_y")["index_x"].values

true_profs = true_profs[subset_inds]
pred_profs = pred_profs[subset_inds]
pred_coords = all_pred_coords[subset_inds]

# Make sure the coordinates all match
assert np.all(pred_coords == shap_coords)

# Import the TF-MoDISco results object
tfm_obj = import_tfmodisco_results(tfm_results_path, hyp_scores, one_hot_seqs, shap_score_center_size)

Plot some SHAP score tracks

Plot the central region of some randomly selected actual importance scores

plot_slice = slice(int(shap_score_center_size / 4), int(3 * shap_score_center_size / 4))
for index in np.random.choice(hyp_scores.shape[0], size=5, replace=False):
    viz_sequence.plot_weights((hyp_scores[index] * one_hot_seqs[index])[plot_slice], subticks_frequency=100)

Plot TF-MoDISco results

Plot all motifs by metacluster

motif_pfms, motif_hcwms, motif_cwms = [], [], []  # Save the trimmed PFMs, hCWMs, and CWMs
motif_pfms_short = []  # PFMs that are even more trimmed (for TOMTOM)
num_seqlets = []  # Number of seqlets for each motif
motif_seqlets = []  # Save seqlets of each motif
if tfm_motifs_cache_dir:
    motif_hdf5 = h5py.File(os.path.join(tfm_motifs_cache_dir, "all_motifs.h5"), "w")
metaclusters = tfm_obj.metacluster_idx_to_submetacluster_results
num_metaclusters = len(metaclusters.keys())
for metacluster_i, metacluster_key in enumerate(metaclusters.keys()):
    metacluster = metaclusters[metacluster_key]
    display(vdomh.h3("Metacluster %d/%d" % (metacluster_i + 1, num_metaclusters)))
    patterns = metacluster.seqlets_to_patterns_result.patterns
    if not patterns:
        break
    motif_pfms.append([])
    motif_hcwms.append([])
    motif_cwms.append([])
    motif_pfms_short.append([])
    num_seqlets.append([])
    motif_seqlets.append([])
    num_patterns = len(patterns)
    for pattern_i, pattern in enumerate(patterns):
        seqlets = pattern.seqlets
        display(vdomh.h4("Pattern %d/%d" % (pattern_i + 1, num_patterns)))
        display(vdomh.p("%d seqlets" % len(seqlets)))
        
        pfm = pattern["sequence"].fwd
        hcwm = pattern["task0_hypothetical_contribs"].fwd
        cwm = pattern["task0_contrib_scores"].fwd
        
        pfm_fig = viz_sequence.plot_weights(pfm, subticks_frequency=10, return_fig=True)
        hcwm_fig = viz_sequence.plot_weights(hcwm, subticks_frequency=10, return_fig=True)
        cwm_fig = viz_sequence.plot_weights(cwm, subticks_frequency=10, return_fig=True)
        pfm_fig.tight_layout()
        hcwm_fig.tight_layout()
        cwm_fig.tight_layout()
        
        motif_table = vdomh.table(
            vdomh.tr(
                vdomh.td("Sequence (PFM)"),
                vdomh.td(figure_to_vdom_image(pfm_fig))
            ),
            vdomh.tr(
                vdomh.td("Hypothetical contributions (hCWM)"),
                vdomh.td(figure_to_vdom_image(hcwm_fig))
            ),
            vdomh.tr(
                vdomh.td("Actual contributions (CWM)"),
                vdomh.td(figure_to_vdom_image(cwm_fig))
            )
        )
        display(motif_table)
        plt.close("all")  # Remove all standing figures
        
        # Trim motif based on information content
        short_trimmed_pfm = trim_motif_by_ic(pfm, pfm)
        motif_pfms_short[-1].append(short_trimmed_pfm)
        
        # Expand trimming to +/- 4bp on either side
        trimmed_pfm = trim_motif_by_ic(pfm, pfm, pad=4)
        trimmed_hcwm = trim_motif_by_ic(pfm, hcwm, pad=4)
        trimmed_cwm = trim_motif_by_ic(pfm, cwm, pad=4)
        
        motif_pfms[-1].append(trimmed_pfm)
        motif_hcwms[-1].append(trimmed_hcwm)
        motif_cwms[-1].append(trimmed_cwm)
        
        num_seqlets[-1].append(len(seqlets))
        
        seqlet_true_profs, seqlet_pred_profs, seqlet_seqs, seqlet_hyps, seqlet_coords = extract_profiles_and_coords(
            seqlets, one_hot_seqs, hyp_scores, true_profs, pred_profs, pred_coords,
            input_length, profile_length, shap_score_center_size,
            profile_display_center_size, task_index=task_index
        )
        
        motif_seqlets[-1].append((seqlet_seqs, seqlet_hyps))

        assert np.allclose(np.sum(seqlet_seqs, axis=0) / len(seqlet_seqs), pattern["sequence"].fwd)
        # ^Sanity check: PFM derived from seqlets match the PFM stored in the pattern
        prof_fig = plot_profiles(
            # Flatten to NT x O x 2
            np.reshape(seqlet_true_profs, (-1, seqlet_true_profs.shape[2], seqlet_true_profs.shape[3])),
            np.reshape(seqlet_pred_profs, (-1, seqlet_pred_profs.shape[2], seqlet_pred_profs.shape[3]))
        )
        
        summit_dists = get_summit_distances(seqlet_coords, peak_table)
        dist_fig = plot_summit_dists(summit_dists)
        
        if tfm_motifs_cache_dir:
            # Save results and figures
            motif_id = "%d_%d" % (metacluster_i, pattern_i)
            np.savez_compressed(
                os.path.join(tfm_motifs_cache_dir, motif_id + "_seqlets"),
                seqlet_true_profs=seqlet_true_profs, seqlet_pred_profs=seqlet_pred_profs,
                seqlet_seqs=seqlet_seqs, seqlet_hyps=seqlet_hyps, seqlet_coords=seqlet_coords,
                summit_dists=summit_dists
            )
            prof_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_profiles.png"))
            dist_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_summitdists.png"))
            pfm_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_pfm_full.png"))
            hcwm_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_hcwm_full.png"))
            cwm_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_cwm_full.png"))
            motif_dset = motif_hdf5.create_group(motif_id)
            motif_dset.create_dataset("pfm_full", data=pfm, compression="gzip")
            motif_dset.create_dataset("hcwm_full", data=hcwm, compression="gzip")
            motif_dset.create_dataset("cwm_full", data=cwm, compression="gzip")
            motif_dset.create_dataset("pfm_trimmed", data=trimmed_pfm, compression="gzip")
            motif_dset.create_dataset("hcwm_trimmed", data=trimmed_hcwm, compression="gzip")
            motif_dset.create_dataset("cwm_trimmed", data=trimmed_cwm, compression="gzip")
            motif_dset.create_dataset("pfm_short_trimmed", data=short_trimmed_pfm, compression="gzip")
if tfm_motifs_cache_dir:
    motif_hdf5.close()

Metacluster 1/1

Pattern 1/7

5646 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Pattern 2/7

3908 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Pattern 3/7

3431 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Pattern 4/7

2508 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Pattern 5/7

129 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Pattern 6/7

127 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Pattern 7/7

52 seqlets

Sequence (PFM)
Hypothetical contributions (hCWM)
Actual contributions (CWM)

Summary of motifs

Motifs are trimmed based on information content, and presented in descending order by number of supporting seqlets. The motifs are separated by metacluster. The motifs are presented as hCWMs. The forward orientation is defined as the orientation that is richer in purines.

colgroup = vdomh.colgroup(
    vdomh.col(style={"width": "5%"}),
    vdomh.col(style={"width": "5%"}),
    vdomh.col(style={"width": "45%"}),
    vdomh.col(style={"width": "45%"})
)
header = vdomh.thead(
    vdomh.tr(
        vdomh.th("#", style={"text-align": "center"}),
        vdomh.th("Seqlets", style={"text-align": "center"}),
        vdomh.th("Forward", style={"text-align": "center"}),
        vdomh.th("Reverse", style={"text-align": "center"})
    )
)

for i in range(len(motif_hcwms)):
    display(vdomh.h3("Metacluster %d/%d" % (i + 1, num_metaclusters)))
    body = []
    for j in range(len(motif_hcwms[i])):
        motif = motif_hcwms[i][j]
        if np.sum(motif[:, [0, 2]]) > 0.5 * np.sum(motif):
            # Forward is purine-rich, reverse-complement is pyrimidine-rich
            f, rc = motif, np.flip(motif, axis=(0, 1))
        else:
            f, rc = np.flip(motif, axis=(0, 1)), motif
            
        f_fig = viz_sequence.plot_weights(f, figsize=(20, 4), return_fig=True)
        f_fig.tight_layout()
        rc_fig = viz_sequence.plot_weights(rc, figsize=(20, 4), return_fig=True)
        rc_fig.tight_layout()
        
        if tfm_motifs_cache_dir:
            # Save results and figures
            motif_id = "%d_%d" % (i, j)
            f_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_hcwm_trimmed_fwd.png"))
            rc_fig.savefig(os.path.join(tfm_motifs_cache_dir, motif_id + "_hcwm_trimmed_rev.png"))

        body.append(
            vdomh.tr(
                vdomh.td(str(j + 1)),
                vdomh.td(str(num_seqlets[i][j])),
                vdomh.td(figure_to_vdom_image(f_fig)),
                vdomh.td(figure_to_vdom_image(rc_fig))
            )
        )
    display(vdomh.table(colgroup, header, vdomh.tbody(*body)))
    plt.close("all")

Metacluster 1/1

#	Seqlets	Forward	Reverse
1	5646
2	3908
3	3431
4	2508
5	129
6	127
7	52

Top TOMTOM matches for each motif

Here, the TF-MoDISco motifs are plotted as hCWMs, but the TOMTOM matches are shown as PWMs.

num_matches_to_keep = 10
num_matches_to_show = 5

header = vdomh.thead(
    vdomh.tr(
        vdomh.th("Motif ID", style={"text-align": "center"}),
        vdomh.th("q-val", style={"text-align": "center"}),
        vdomh.th("PWM", style={"text-align": "center"})
    )
)

for i in range(len(motif_pfms)):
    display(vdomh.h3("Metacluster %d/%d" % (i + 1, num_metaclusters)))
    
    # Compute TOMTOM matches for all motifs in the metacluster at once
    out_dir = os.path.join(tfm_motifs_cache_dir, "tomtom", "metacluster_%d" % i) if tfm_motifs_cache_dir else None
    tomtom_matches = match_motifs_to_database(motif_pfms_short[i], top_k=num_matches_to_keep, temp_dir=out_dir)
    
    for j in range(len(motif_pfms[i])):
        display(vdomh.h4("Motif %d/%d" % (j + 1, len(motif_pfms[i]))))
        viz_sequence.plot_weights(motif_hcwms[i][j])
    
        body = []
        for k, (match_name, match_pfm, match_qval) in enumerate(tomtom_matches[j]):
            fig = viz_sequence.plot_weights(pfm_to_pwm(match_pfm), return_fig=True)
            fig.tight_layout()
            if k < num_matches_to_show:
                body.append(
                    vdomh.tr(
                        vdomh.td(match_name),
                        vdomh.td(str(match_qval)),
                        vdomh.td(figure_to_vdom_image(fig))
                    )
                )
                if tfm_motifs_cache_dir:
                    # Save results and figures
                    motif_id = "%d_%d" % (i, j)
                    fig.savefig(os.path.join(out_dir, motif_id + ("_hit-%d.png" % (k + 1))))
            else:
                body.append(
                    vdomh.tr(
                        vdomh.td(match_name),
                        vdomh.td(str(match_qval)),
                        vdomh.td("Not shown")
                    )
                )
        if not body:
            display(vdomh.p("No TOMTOM matches passing threshold"))
        else:
            display(vdomh.table(header, vdomh.tbody(*body)))
        plt.close("all")

Metacluster 1/1

Motif 1/7

Motif ID	q-val	PWM
MA0104.4_MYCN	0.000234888
MA0147.3_MYC	0.000234888
MYCN_HUMAN.H11MO.0.A	0.000234888
MAX_HUMAN.H11MO.0.A	0.0017433
MA0059.1_MAX::MYC	0.00181594
MXI1_HUMAN.H11MO.0.A	0.00280173	Not shown
CR3L1_HUMAN.H11MO.0.D	0.00355614	Not shown
MA0004.1_Arnt	0.0045587	Not shown
MXI1_HUMAN.H11MO.1.A	0.00496983	Not shown
MA0608.1_Creb3l2	0.00655202	Not shown

Motif 2/7

Motif ID	q-val	PWM
MA1099.2_HES1	0.14765899999999998
MYC_HUMAN.H11MO.0.A	0.14765899999999998
MA0746.2_SP3	0.14765899999999998
MA0608.1_Creb3l2	0.14765899999999998
MA0626.1_Npas2	0.14765899999999998
MA0079.4_SP1	0.14765899999999998	Not shown
MYCN_HUMAN.H11MO.0.A	0.14765899999999998	Not shown
MA0747.1_SP8	0.14765899999999998	Not shown
BMAL1_HUMAN.H11MO.0.A	0.14765899999999998	Not shown
MA1513.1_KLF15	0.14765899999999998	Not shown

Motif 3/7

Motif ID	q-val	PWM
MA1108.2_MXI1	0.0016257
MYC_HUMAN.H11MO.0.A	0.0193988
MA0147.3_MYC	0.0193988
MA0668.1_NEUROD2	0.0193988
MA0104.4_MYCN	0.0193988
MA0058.3_MAX	0.0193988	Not shown
MA0059.1_MAX::MYC	0.0227547	Not shown
MA1568.1_TCF21(var.2)	0.029812400000000003	Not shown
MA0825.1_MNT	0.029812400000000003	Not shown
MAX_HUMAN.H11MO.0.A	0.032626999999999996	Not shown

Motif 4/7

Motif ID	q-val	PWM
CTCF_HUMAN.H11MO.0.A	3.49485e-16
MA0139.1_CTCF	4.71717e-13
CTCFL_HUMAN.H11MO.0.A	2.84301e-08
MA1102.2_CTCFL	3.2653899999999996e-05
SNAI1_HUMAN.H11MO.0.C	0.15193099999999998
MA1568.1_TCF21(var.2)	0.15193099999999998	Not shown
MA1638.1_HAND2	0.20521399999999998	Not shown
MA0155.1_INSM1	0.253689	Not shown
KLF8_HUMAN.H11MO.0.C	0.285926	Not shown
ZIC3_HUMAN.H11MO.0.B	0.321162	Not shown

Motif 5/7

Motif ID	q-val	PWM
SP2_HUMAN.H11MO.0.A	0.00862481
SP3_HUMAN.H11MO.0.B	0.00862481
MA1650.1_ZBTB14	0.0411422
SP1_HUMAN.H11MO.0.A	0.0411422
MXI1_HUMAN.H11MO.0.A	0.0411422
SP1_HUMAN.H11MO.1.A	0.05604069999999999	Not shown
KLF3_HUMAN.H11MO.0.B	0.08319689999999999	Not shown
PATZ1_HUMAN.H11MO.0.C	0.08375650000000001	Not shown
CTCFL_HUMAN.H11MO.0.A	0.08375650000000001	Not shown
KLF9_HUMAN.H11MO.0.C	0.08375650000000001	Not shown

Motif 6/7

Motif ID	q-val	PWM
NRF1_HUMAN.H11MO.0.A	0.131946
HES7_HUMAN.H11MO.0.D	0.131946
MA1099.2_HES1	0.131946
MA0506.1_NRF1	0.131946
MA0821.1_HES5	0.131946
HES5_HUMAN.H11MO.0.D	0.131946	Not shown
MA0616.2_HES2	0.131946	Not shown
MA0822.1_HES7	0.136925	Not shown
ZBT14_HUMAN.H11MO.0.C	0.136925	Not shown
AP2D_HUMAN.H11MO.0.D	0.136925	Not shown

Motif 7/7

Motif ID	q-val	PWM
MA0076.2_ELK4	0.00690137
WT1_HUMAN.H11MO.0.C	0.00690137
ZN467_HUMAN.H11MO.0.C	0.00690137
SP3_HUMAN.H11MO.0.B	0.00690137
KLF3_HUMAN.H11MO.0.B	0.00743804
PATZ1_HUMAN.H11MO.0.C	0.0106088	Not shown
MAZ_HUMAN.H11MO.0.A	0.0115549	Not shown
ETV1_HUMAN.H11MO.0.A	0.0115549	Not shown
FEV_HUMAN.H11MO.0.B	0.0115549	Not shown
ELK4_HUMAN.H11MO.0.A	0.0127195	Not shown

Sample of seqlets supporting each motif

Here, the motifs are presented as hCWMs, along with the actual importance scores of a random sample of seqlets that support the motif.

num_seqlets_to_show = 10

colgroup = vdomh.colgroup(
    vdomh.col(style={"width": "50%"}),
    vdomh.col(style={"width": "50%"})
)

header = vdomh.thead(
    vdomh.tr(
        vdomh.th("Motif hCWM", style={"text-align": "center"}),
        vdomh.th("Seqlets", style={"text-align": "center"})
    )
)

for i in range(len(motif_hcwms)):
    display(vdomh.h3("Metacluster %d/%d" % (i + 1, num_metaclusters)))
    
    for j in range(len(motif_hcwms[i])):
        display(vdomh.h4("Motif %d/%d" % (j + 1, len(motif_hcwms[i]))))
        
        motif_fig = viz_sequence.plot_weights(motif_hcwms[i][j], figsize=(20, 4), return_fig=True)
        motif_fig.tight_layout()
        
        seqlet_seqs, seqlet_hyps = motif_seqlets[i][j]
        
        sample_size = min(num_seqlets_to_show, len(seqlet_seqs))
        sample_inds = np.random.choice(len(seqlet_seqs), size=sample_size, replace=False)
        sample = []
        for k in sample_inds:
            fig = viz_sequence.plot_weights(seqlet_hyps[k] * seqlet_seqs[k], subticks_frequency=10, return_fig=True)
            fig.tight_layout()
            sample.append(figure_to_vdom_image(fig))
        body = vdomh.tbody(vdomh.tr(vdomh.td(figure_to_vdom_image(motif_fig)), vdomh.td(*sample)))
        display(vdomh.table(colgroup, header, body))
        plt.close("all")

Metacluster 1/1

Motif 1/7

Motif hCWM	Seqlets

Motif 2/7

Motif hCWM	Seqlets

Motif 3/7

Motif hCWM	Seqlets

Motif 4/7

Motif hCWM	Seqlets

Motif 5/7

Motif hCWM	Seqlets

Motif 6/7

Motif hCWM	Seqlets

Motif 7/7

Motif hCWM	Seqlets