lasseignelab · lizzyjoan · May 1, 2026 · Nov 24, 2025 · Nov 24, 2025 · Nov 24, 2025
diff --git a/.gitignore b/.gitignore
@@ -241,6 +241,10 @@ bin/conda/work/
 bin/conda/work/*
 bin/conda/.nextflow/
 bin/container/*
+<<<<<<< HEAD
+bin/docker/
+src/lint_and_style.sh
+src/lint_and_style.sh
 results/permutation_analysis/aggregate_gene_exp_mat.csv
 Rplots.pdf
 src/lint_and_style.sh

diff --git a/config/environments/lasseignelab.sh b/config/environments/lasseignelab.sh
@@ -12,6 +12,9 @@ cap_data_link "$LASSEIGNE_LAB_PATH/PROJECT_dir/setbp1_hd/tanaka"
 # Wong et al 2025 source data to validate DEGs
 cap_data_link "$LASSEIGNE_LAB_PATH/PROJECT_dir/setbp1_hd/wong"
 
+# DailyMed drug safety data
+cap_data_link "$LASSEIGNE_LAB_PATH/PROJECT_dir/setbp1_hd/dailymed" 
+
 # Brainspan and GTEx counts data
 cap_data_link "$LASSEIGNE_LAB_PATH/DATASET_dir/brainspan"
 cap_data_link "$LASSEIGNE_LAB_PATH/DATASET_dir/gtex"
diff --git a/results/figures/approved_drugs_braincells_full.png b/results/figures/approved_drugs_braincells_full.png
diff --git a/results/figures/approved_drugs_braincells_multiplecells.png b/results/figures/approved_drugs_braincells_multiplecells.png
diff --git a/results/figures/approved_drugs_braincells_topWTCS.png b/results/figures/approved_drugs_braincells_topWTCS.png
diff --git a/results/figures/approved_drugs_top_moa.png b/results/figures/approved_drugs_top_moa.png
diff --git a/results/figures/drug_targets_histogram.png b/results/figures/drug_targets_histogram.png
diff --git a/results/figures/figure_5.png b/results/figures/figure_5.png
diff --git a/results/figures/full_drug_safety_anno_table.csv b/results/figures/full_drug_safety_anno_table.csv
diff --git a/results/figures/setbp1_converge_drug_anno_table.csv b/results/figures/setbp1_converge_drug_anno_table.csv
diff --git a/results/figures/setbp1_drug_shared_targets_trialphase.png b/results/figures/setbp1_drug_shared_targets_trialphase.png
diff --git a/results/figures/setbp1_drugs_shared.png b/results/figures/setbp1_drugs_shared.png
diff --git a/results/figures/setbp1_drugs_shared_phase.png b/results/figures/setbp1_drugs_shared_phase.png
diff --git a/results/figures/setbp1_drugs_shared_targets.png b/results/figures/setbp1_drugs_shared_targets.png
diff --git a/results/signature_reversion/approved_drugs_res.csv b/results/signature_reversion/approved_drugs_res.csv
diff --git a/results/signature_reversion/approveddrugs_target_exp_heatmap.png b/results/signature_reversion/approveddrugs_target_exp_heatmap.png
diff --git a/results/signature_reversion/brain_drugs_filtered.csv b/results/signature_reversion/brain_drugs_filtered.csv
diff --git a/results/signature_reversion/drug_target_degs.csv b/results/signature_reversion/drug_target_degs.csv
diff --git a/results/signature_reversion/lincs_anno_braincell_results.rds b/results/signature_reversion/lincs_anno_braincell_results.rds
diff --git a/results/signature_reversion/lincs_braincells_results.csv b/results/signature_reversion/lincs_braincells_results.csv
diff --git a/results/signature_reversion/lincs_results.rds b/results/signature_reversion/lincs_results.rds
diff --git a/results/signature_reversion/setbp1_drug_targets.csv b/results/signature_reversion/setbp1_drug_targets.csv
diff --git a/results/signature_reversion/target_exp_contrasts.png b/results/signature_reversion/target_exp_contrasts.png
diff --git a/results/signature_reversion/target_exp_heatmap.png b/results/signature_reversion/target_exp_heatmap.png
diff --git a/src/11_job_sigsearch.sh b/src/11_job_sigsearch.sh
@@ -0,0 +1,16 @@
+#!/bin/bash
+
+#SBATCH --ntasks=1
+#SBATCH --mem-per-cpu=32G
+#SBATCH --partition=express
+
+cap_container -c singularity "lizzyr/sigsearch:1.0.0"
+
+singularity exec --cleanenv \
+                 --containall \
+                 -B "${CAP_PROJECT_PATH}" \
+                 -B "${LASSEIGNE_LAB_PATH}" \
+                 "${CAP_CONTAINER_PATH}"/sigsearch_1.0.0.sif \
+                  Rscript --vanilla "${CAP_PROJECT_PATH}"/src/11_signature_reversion.R \
+                  -i "${CAP_DATA_PATH}"/ \
+                  -o "${CAP_RESULTS_PATH}"/
diff --git a/src/11_signature_reversion.R b/src/11_signature_reversion.R
@@ -0,0 +1,202 @@
+suppressPackageStartupMessages({
+  library(argparse)
+  library(tidyverse)
+  library(recount3)
+  library(signatureSearch)
+  library(ExperimentHub)
+  library(rhdf5)
+})
+
+# create parser object
+parser <- ArgumentParser()
+
+parser$add_argument("-i", "--input_dir",
+  type = "character",
+  metavar = "PATH",
+  help = "Path to data directory"
+)
+parser$add_argument("-o", "--output_dir",
+  type = "character",
+  metavar = "PATH",
+  help = "Path to the results directory to read consensus signature and write signatureSearch results" # nolint
+)
+
+# get command line options, if help option encountered print help and exit
+args <- parser$parse_args()
+
+# source functions file if supplied
+if (!is.null(args$functions)) source(args$functions)
+
+# start timer
+ptm <- proc.time()
+
+#### Set variables for dataset ####
+
+datapath <- args$input_dir
+resdir <- args$output_dir
+filepath <- paste0(resdir, "signature_reversion/")
+querypath <- paste0(resdir, "consensus_signature/")
+cachedir <- paste0(datapath, "ExperimentHubCache")
+message("Location of cache directory: ", cachedir)
+
+# create dirs if needed
+if (!dir.exists(filepath)) dir.create(filepath, recursive = TRUE)
+if (!dir.exists(cachedir)) dir.create(cachedir, recursive = TRUE)
+
+###### EXPERIMENT HUB SET UP ######
+
+setExperimentHubOption("ASK", FALSE)
+setExperimentHubOption("CACHE", cachedir)
+
+# Force a clean load of the hub
+Sys.setenv(EXPERIMENT_HUB_CACHE = cachedir)
+eh <- ExperimentHub(localHub = FALSE)
+
+# lincs 2020 filepath and db exploration
+message("Loading in lincs2 EH db...")
+lincs2 <- eh[["EH7297"]]
+
+# load LINCS perturbation annotations
+data(lincs_pert_info2)
+message("Loading in lincs2 annotation EH db...")
+
+###### READ IN DATA ######
+
+# RRA signature
+rra_degs <- read_csv(
+  paste0(querypath, "rankaggregate_consensus_signature.csv")
+)
+# SETBP1 regulatory targets
+# Curated targets by Sasha Taluri
+setbp1_targets_all <- read_csv(paste0(
+  datapath, "260403_setbp1_targets.csv"
+))
+# Prioritized targets from coexpression & pathway analysis
+setbp1_targets_top <- read_csv(paste0(
+  resdir, "permutation_analysis/prioritized_setbp1_targets.csv"
+))
+
+###### PREPARE QUERY SIGNATURE ######
+
+# genes available in LINCS DB
+db_genes <- h5read(lincs2, "rownames", drop = TRUE)
+
+# up vs down DEGs in signature
+degs_up <- rra_degs %>%
+  mutate(entrez = as.character(entrez)) %>%
+  filter(meta_logFC > 0.1) %>%
+  filter(entrez %in% db_genes) %>%
+  slice_min(Score, n = 100)
+cat("Up genes summary... ", summary(degs_up), "\n", sep = " ")
+
+degs_down <- rra_degs %>%
+  mutate(entrez = as.character(entrez)) %>%
+  filter(meta_logFC < -0.1) %>%
+  filter(entrez %in% db_genes) %>%
+  slice_min(Score, n = 100)
+cat("Down genes summary... ", summary(degs_down), "\n", sep = " ")
+
+###### SIGNATURE REVERSION ######
+
+# query LINCS signatures
+qsig_lincs <- qSig(
+  query = list(upset = degs_up$entrez, downset = degs_down$entrez),
+  gess_method = "LINCS",
+  refdb = lincs2
+)
+
+# calculate signature reversion metrics of query vs LINCS DB
+lincs <- gess_lincs(qsig_lincs, sortby = "WTCS", tau = FALSE, workers = 1)
+str(lincs)
+write_rds(lincs, paste0(filepath, "lincs_results.rds"))
+
+# filter for brain/CNS-derived/relevant cells (HEK cells often used in neuro)
+brain_cells <- c("NEU", "SHSY5Y", "NPC", "HEK293", "HEK293T")
+lincs_brain <- lincs@result %>%
+  filter(cell %in% brain_cells) %>%
+  filter(WTCS < 0) %>%
+  arrange(WTCS)
+head(lincs_brain)
+write_csv(lincs_brain, paste0(filepath, "lincs_braincells_results.csv"))
+
+# annotate drugs
+brain_drugs_anno <- lincs_brain %>%
+  left_join(lincs_pert_info2, by = c("pert" = "pert_id"))
+cat(
+  "Number of unique reversion drugs in brain-relevant cells: ",
+  length(unique(brain_drugs_anno$pert)), "\n",
+  sep = " "
+)
+
+write_rds(
+  brain_drugs_anno, paste0(filepath, "lincs_anno_braincell_results.rds")
+)
+
+###### TARGET ANALYSIS ######
+# separate targets by rows to compare drug targets
+brain_drug_targets <- brain_drugs_anno %>%
+  dplyr::select(
+    pert, pref_name, max_phase, cell, WTCS,
+    WTCS_Pval, mergeTargets, Target_pathway
+  ) %>%
+  separate_longer_delim(mergeTargets, delim = "; ") %>%
+  filter(!is.na(mergeTargets)) %>%
+  add_count(mergeTargets) %>%
+  # annotate targets that are prioritized SETBP1 regulatory targets
+  mutate(
+    SETBP1_target = ifelse(
+      mergeTargets %in% setbp1_targets_all$Target, TRUE, FALSE
+    ),
+    SETBP1_top_target = ifelse(
+      mergeTargets %in% setbp1_targets_top$Target_Gene, TRUE, FALSE
+    )
+  )
+
+# filter for drugs that target any SETBP1 regulatory targets
+drug_targets_setbp1 <- brain_drug_targets %>%
+  filter(SETBP1_target == TRUE) %>%
+  distinct()
+str(drug_targets_setbp1)
+print(n = 21, arrange(drug_targets_setbp1, pref_name))
+
+cat(
+  "Drugs that target SETBP1 targets: ",
+  unique(drug_targets_setbp1$pref_name),
+  "\n",
+  sep = " "
+)
+cat(
+  "Drug target SETBP1 targets: ",
+  unique(drug_targets_setbp1$mergeTargets),
+  "\n",
+  sep = " "
+)
+write_csv(
+  drug_targets_setbp1, paste0(filepath, "setbp1_drug_targets.csv")
+)
+
+# full drug annotations for setbp1-target-targeting drugs
+
+setbp1_drugs <- brain_drugs_anno %>%
+  filter(pert %in% drug_targets_setbp1$pert)
+
+
+#### APPROVED/POST-MARKET DRUGS ####
+approved_drugs <- brain_drugs_anno %>%
+  # max WTCS to help with plot ordering
+  group_by(pref_name) %>%
+  mutate(drug_top_WTCS = min(WTCS)) %>%
+  ungroup() %>%
+  add_count(pert) %>%
+  # filter for approved drugs
+  filter(max_phase == 4)
+
+write_csv(approved_drugs, paste0(filepath, "approved_drugs_res.csv"))
+
+
+# end timer
+fptm <- proc.time()
+(fptm[3] / 60)
+
+#  session info
+print(sessionInfo())