Skip to content

Commit

Permalink
Merge branch 'feature/wilms-tumor-06-azimuth' into sjspielman/wilms-0…
Browse files Browse the repository at this point in the history 
…6_update-infercnv
  • Loading branch information
@sjspielman
sjspielman authored Nov 12, 2024
2 parents 93636b1 + c1d0bcc commit 8f15a2d
Show file tree
Hide file tree
Showing 8 changed files with 219 additions and 284 deletions.
7 changes: 6 additions & 1 deletion analyses/cell-type-wilms-tumor-06/00_run_workflow.sh
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,9 @@
# steps that do not directly contribute to final cell type annotations
# should be run. Setting RUN_EXPLORATORY=1 will run those steps.
# This variable is 0 by default.
# 3. The THREADS variable controls how many cores are used for inference
# with copyKAT, which is an exploratory step in the workflow.
# The variable is 32 by default.
#
# Default usage:
# ./00_run_workflow.sh
Expand All @@ -24,6 +27,7 @@ set -euo pipefail

IS_CI=${TESTING:-0}
RUN_EXPLORATORY=${RUN_EXPLORATORY:-0}
THREADS=${THREADS:-32}
project_id="SCPCP000006"

# Ensure script is being run from its directory
Expand Down Expand Up @@ -130,7 +134,8 @@ if [[ $RUN_EXPLORATORY -eq 1 ]]; then

# Run infercnv and copykat for a selection of samples
# This script calls scripts/05_copyKAT.R and scripts/06_infercnv.R
Rscript scripts/explore-cnv-methods.R
# By default, copyKAT as called by this script uses 32 cores
THREADS=${THREADS} TESTING=${IS_CI} ./scripts/explore-cnv-methods.sh

fi

Expand Down
89 changes: 51 additions & 38 deletions analyses/cell-type-wilms-tumor-06/notebook_template/05_copykat_exploration.Rmd
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@ title: "Copykat CNV results exploration for `r params$sample_id`"
author: "Maud PLASCHKA"
date: "`r Sys.Date()`"
params:
sample_id: "SCPCS000194"
sample_id: "SCPCS000179"
seed: 12345
output:
html_document:
Expand All @@ -16,27 +16,24 @@ output:


```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE,
message=FALSE,
warnings=FALSE
knitr::opts_chunk$set(
echo = TRUE,
message = FALSE,
warnings = FALSE
)
```


## Introduction


```{r}
subdiagnosis <- readr::read_tsv(
file.path("..", "..", "..", "data", "current", "SCPCP000006", "single_cell_metadata.tsv"),
show_col_types = FALSE
) |>
dplyr::filter(scpca_sample_id == params$sample_id) |>
dplyr::pull(subdiagnosis)
```

This notebook explores using [`CopyKAT`](https://github.com/navinlabcode/copykat) to estimate tumor and normal cells in `r params$sample_id` from SCPCP000006.
Expand All @@ -59,7 +56,7 @@ These results are read into this notebook and used to:

```{r packages, message=FALSE, warning=FALSE}
library(Seurat)
library(SCpubr) # for plotting
library(SCpubr) # for plotting
library(tidyverse)
library(patchwork)
Expand Down Expand Up @@ -87,14 +84,33 @@ In this notebook, we are working with the Wilms tumor sample defined in `r param
We work with the pre-processed and labeled `Seurat` object that is the output of `02b_label-transfer_fetal_kidney_reference_Stewart.Rmd` saved in the `results` directory.

```{r}
result_dir <- file.path(module_base, "results", params$sample_id)
predictions_paths <- list()
full_ck_result_paths <- list()
for(ref_value in c("ref", "noref")){
for(distance_value in c("euclidean", "spearman")){
predictions_file <- glue::glue("05_copykat_", {params$sample_id}, "_",ref_value,"_distance-", distance_value, "_copykat_prediction.txt")
full_ck_result_file <- glue::glue("05_copykat_", {params$sample_id}, "_",ref_value,"_distance-", distance_value, "_copykat_CNA_results.txt")
for (ref_value in c("ref", "noref")) {
for (distance_value in c("euclidean", "spearman")) {
predictions_file <- glue::glue(
"05_copykat_",
{
params$sample_id
},
"_",
ref_value,
"_distance-",
distance_value,
"_copykat_prediction.txt"
)
full_ck_result_file <- glue::glue(
"05_copykat_",
{
params$sample_id
},
"_",
ref_value,
"_distance-",
distance_value,
"_copykat_CNA_results.txt"
)
predictions_paths[[glue::glue(ref_value, "_", distance_value)]] <- file.path(result_dir, "05_copyKAT", ref_value, distance_value, predictions_file)
full_ck_result_paths[[glue::glue(ref_value, "_", distance_value)]] <- file.path(result_dir, "05_copyKAT", ref_value, distance_value, full_ck_result_file)
Expand All @@ -121,7 +137,7 @@ Here we defined function that will be used multiple time all along the notebook.
### Load the pre-processed `Seurat` object
```{r load, message=FALSE, warning=FALSE}
# open the processed rds object
srat <- readRDS(file.path(result_dir, paste0("02b-fetal_kidney_label-transfer_", params$sample_id,".Rds")))
srat <- readRDS(file.path(result_dir, paste0("02b-fetal_kidney_label-transfer_", params$sample_id, ".Rds")))
DefaultAssay(srat) <- "SCT"
```

Expand All @@ -134,21 +150,21 @@ Below we look at the heatmaps produced by `CopyKAT`.

##### Euclidean distance

![](`r file.path(result_dir, "05_copyKAT", "noref", "euclidean", glue::glue("05_copykat_", params$sample_id, "_noref_distance-euclidean_copykat_heatmap.png"))`)
![](`r file.path(result_dir, "05_copyKAT", "noref", "euclidean", glue::glue("05_copykat_", params$sample_id, "_noref_distance-euclidean_copykat_heatmap.jpg"))`)

##### Spearman distance

![](`r file.path(result_dir, "05_copyKAT", "noref", "spearman", glue::glue("05_copykat_", params$sample_id, "_noref_distance-spearman_copykat_heatmap.png"))`)
![](`r file.path(result_dir, "05_copyKAT", "noref", "spearman", glue::glue("05_copykat_", params$sample_id, "_noref_distance-spearman_copykat_heatmap.jpg"))`)

#### Heatmap with endothelial cells as reference

##### Euclidean distance

![](`r file.path(result_dir, "05_copyKAT", "ref", "euclidean", glue::glue("05_copykat_", params$sample_id, "_ref_distance-euclidean_copykat_heatmap.png"))`)
![](`r file.path(result_dir, "05_copyKAT", "ref", "euclidean", glue::glue("05_copykat_", params$sample_id, "_ref_distance-euclidean_copykat_heatmap.jpg"))`)

##### Spearman distance

![](`r file.path(result_dir, "05_copyKAT", "ref", "spearman", glue::glue("05_copykat_", params$sample_id, "_ref_distance-spearman_copykat_heatmap.png"))`)
![](`r file.path(result_dir, "05_copyKAT", "ref", "spearman", glue::glue("05_copykat_", params$sample_id, "_ref_distance-spearman_copykat_heatmap.jpg"))`)


#### UMAP
Expand All @@ -158,19 +174,17 @@ We show a side by side UMAP with results from running `CopyKAT` both with and wi

```{r}
# read in ck predictions from both reference types (no_normal and with_normal)
ck_results_df <- predictions_paths |>
purrr::map(readr::read_tsv) |>
ck_results_df <- predictions_paths |>
purrr::map(readr::read_tsv) |>
dplyr::bind_rows(.id = "reference_used")
# get umap coordinate
umap_df <- srat[["umap"]]@cell.embeddings |>
as.data.frame() |>
rownames_to_column('barcodes')
rownames_to_column("barcodes")
cnv_df <- umap_df |>
dplyr::left_join(ck_results_df, by = c("barcodes" = "cell.names"))
dplyr::left_join(ck_results_df, by = c("barcodes" = "cell.names"))
```

```{r}
Expand Down Expand Up @@ -201,27 +215,26 @@ We might expect that tumor cells would show an increased estimated copy number i
```{r}
# read in full gene by cell copy number detection results
full_ck_results_df <- full_ck_result_paths |>
purrr::map(readr::read_tsv) |>
dplyr::bind_rows(.id = "reference_used")
purrr::map(readr::read_tsv) |>
dplyr::bind_rows(.id = "reference_used")
# for every cell, calculate the mean detection level across all genes in a given chromosome
full_cnv_df <- full_ck_results_df |>
tidyr::pivot_longer(
tidyr::pivot_longer(
cols = -c(
reference_used,
chrom
reference_used,
chrom
),
names_to = "barcodes",
values_to = "cnv_detection"
) |>
dplyr::group_by(chrom, barcodes, reference_used) |>
dplyr::summarise(mean_cnv_detection = mean(cnv_detection))
) |>
dplyr::group_by(chrom, barcodes, reference_used) |>
dplyr::summarise(mean_cnv_detection = mean(cnv_detection))
# join with cnv info
cnv_df <- cnv_df |>
dplyr::left_join(full_cnv_df, by = c("barcodes", "reference_used")) |>
dplyr::filter(!is.na(chrom))
dplyr::left_join(full_cnv_df, by = c("barcodes", "reference_used")) |>
dplyr::filter(!is.na(chrom))
```
Let's look at the distribution of CNV estimation in cells that are called aneuploid and diploid by `CopyKAT`.
Expand Down
10 changes: 0 additions & 10 deletions analyses/cell-type-wilms-tumor-06/renv.lock
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3165,16 +3165,6 @@
],
"Hash": "8954069286b4b2b0d023d1b288dce978"
},
"jpeg": {
"Package": "jpeg",
"Version": "0.1-10",
"Source": "Repository",
"Repository": "RSPM",
"Requirements": [
"R"
],
"Hash": "031a0b683d001a7519202f0628fc0358"
},
"jquerylib": {
"Package": "jquerylib",
"Version": "0.1.4",
Expand Down
105 changes: 9 additions & 96 deletions analyses/cell-type-wilms-tumor-06/renv/activate.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
local({

# the requested version of renv
version <- "1.0.8"
version <- "1.0.7"
attr(version, "sha") <- NULL

# the project directory
Expand Down Expand Up @@ -98,66 +98,6 @@ local({
unloadNamespace("renv")

# load bootstrap tools
ansify <- function(text) {
if (renv_ansify_enabled())
renv_ansify_enhanced(text)
else
renv_ansify_default(text)
}

renv_ansify_enabled <- function() {

override <- Sys.getenv("RENV_ANSIFY_ENABLED", unset = NA)
if (!is.na(override))
return(as.logical(override))

pane <- Sys.getenv("RSTUDIO_CHILD_PROCESS_PANE", unset = NA)
if (identical(pane, "build"))
return(FALSE)

testthat <- Sys.getenv("TESTTHAT", unset = "false")
if (tolower(testthat) %in% "true")
return(FALSE)

iderun <- Sys.getenv("R_CLI_HAS_HYPERLINK_IDE_RUN", unset = "false")
if (tolower(iderun) %in% "false")
return(FALSE)

TRUE

}

renv_ansify_default <- function(text) {
text
}

renv_ansify_enhanced <- function(text) {

# R help links
pattern <- "`\\?(renv::(?:[^`])+)`"
replacement <- "`\033]8;;ide:help:\\1\a?\\1\033]8;;\a`"
text <- gsub(pattern, replacement, text, perl = TRUE)

# runnable code
pattern <- "`(renv::(?:[^`])+)`"
replacement <- "`\033]8;;ide:run:\\1\a\\1\033]8;;\a`"
text <- gsub(pattern, replacement, text, perl = TRUE)

# return ansified text
text

}

renv_ansify_init <- function() {

envir <- renv_envir_self()
if (renv_ansify_enabled())
assign("ansify", renv_ansify_enhanced, envir = envir)
else
assign("ansify", renv_ansify_default, envir = envir)

}

`%||%` <- function(x, y) {
if (is.null(x)) y else x
}
Expand Down Expand Up @@ -202,10 +142,7 @@ local({
# compute common indent
indent <- regexpr("[^[:space:]]", lines)
common <- min(setdiff(indent, -1L)) - leave
text <- paste(substring(lines, common), collapse = "\n")

# substitute in ANSI links for executable renv code
ansify(text)
paste(substring(lines, common), collapse = "\n")

}

Expand Down Expand Up @@ -369,11 +306,7 @@ local({
)

if ("headers" %in% names(formals(utils::download.file)))
{
headers <- renv_bootstrap_download_custom_headers(url)
if (length(headers) && is.character(headers))
args$headers <- headers
}
args$headers <- renv_bootstrap_download_custom_headers(url)

do.call(utils::download.file, args)

Expand Down Expand Up @@ -452,22 +385,10 @@ local({
for (type in types) {
for (repos in renv_bootstrap_repos()) {

# build arguments for utils::available.packages() call
args <- list(type = type, repos = repos)

# add custom headers if available -- note that
# utils::available.packages() will pass this to download.file()
if ("headers" %in% names(formals(utils::download.file)))
{
headers <- renv_bootstrap_download_custom_headers(url)
if (length(headers) && is.character(headers))
args$headers <- headers
}

# retrieve package database
db <- tryCatch(
as.data.frame(
do.call(utils::available.packages, args),
utils::available.packages(type = type, repos = repos),
stringsAsFactors = FALSE
),
error = identity
Expand Down Expand Up @@ -549,31 +470,23 @@ local({

}

renv_bootstrap_github_token <- function() {
for (envvar in c("GITHUB_TOKEN", "GITHUB_PAT", "GH_TOKEN")) {
envval <- Sys.getenv(envvar, unset = NA)
if (!is.na(envval))
return(envval)
}
}

renv_bootstrap_download_github <- function(version) {

enabled <- Sys.getenv("RENV_BOOTSTRAP_FROM_GITHUB", unset = "TRUE")
if (!identical(enabled, "TRUE"))
return(FALSE)

# prepare download options
token <- renv_bootstrap_github_token()
if (nzchar(Sys.which("curl")) && nzchar(token)) {
pat <- Sys.getenv("GITHUB_PAT")
if (nzchar(Sys.which("curl")) && nzchar(pat)) {
fmt <- "--location --fail --header \"Authorization: token %s\""
extra <- sprintf(fmt, token)
extra <- sprintf(fmt, pat)
saved <- options("download.file.method", "download.file.extra")
options(download.file.method = "curl", download.file.extra = extra)
on.exit(do.call(base::options, saved), add = TRUE)
} else if (nzchar(Sys.which("wget")) && nzchar(token)) {
} else if (nzchar(Sys.which("wget")) && nzchar(pat)) {
fmt <- "--header=\"Authorization: token %s\""
extra <- sprintf(fmt, token)
extra <- sprintf(fmt, pat)
saved <- options("download.file.method", "download.file.extra")
options(download.file.method = "wget", download.file.extra = extra)
on.exit(do.call(base::options, saved), add = TRUE)
Expand Down
Loading

0 comments on commit 8f15a2d

Please sign in to comment.