OSCR

Drug screen and machine learning predict neuroprotective agents in a preclinical human model of childhood dementia.

Code ↔ Paper

23 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 23 matches
  1. [1] § Methods › Single-nuclei data processing ↔ Rscripts/QualityControl_60D.Rmd, lines 247–290 · score 0.96 · generate cell clusters, ambientRNA, FindClusters, FindNeighbors, RunUMAP, SCTransform
  2. [2] § Methods › Single-nuclei data processing ↔ Rscripts/QualityControl_30D.Rmd, lines 416–450 · score 0.96 · generate cell clusters, ambientRNA, FindClusters, FindNeighbors, RunUMAP, SCTransform
  3. [3] § Methods › Cell type classification › scMayoMap cell annotations ↔ Rscripts/Neuronal_Glial_Types_Analysis_60D.Rmd, lines 347–435 · score 0.85 · FindConservedMarkers, cluster_type2, cluster_type1, scMayoMap, cell annotations, database
  4. [4] § Methods › Cell type classification › scMayoMap cell annotations ↔ Rscripts/Neuronal_Glial_Types_Analysis_30D.Rmd, lines 315–392 · score 0.85 · FindConservedMarkers, cluster_type2, cluster_type1, scMayoMap, cell annotations, database
  5. [5] § Methods › Identification of differentially expressed genes ↔ Rscripts/[Sanfilippo - 30D] Differential Expression Analysis.Rmd, lines 207–275 · score 0.80 · min.diff.pct, logfc.threshold, FindMarkers, identify genes, patients
  6. [6] § Methods › Integration with healthy, human cortical biopsies ↔ Rscripts/QualityControl_60D.Rmd, lines 247–290 · score 0.80 · RunPCA, SCT assay, variable features, SCTransform, Principal component, filtered
  7. [7] § Methods › Machine learning feature models ↔ feature models/sanfilipo_cls_machine_learning.ipynb, lines 35–56 · score 0.75 · models hyperparameters, sub model, recursive feature, cross validation, Hyperopt, MPS IIIA
  8. [8] § Methods › Classification of cycling cells ↔ Rscripts/Neuronal_Glial_Types_Analysis_30D.Rmd, lines 200–267 · score 0.73 · cell cycle scores, NonCycling, AddModuleScore, G1S, G2M, gene
  9. [9] § Methods › Cell type classification › UMAP and Seurat cluster-based cell classification ↔ Rscripts/Neuronal_Glial_Types_Analysis_30D.Rmd, lines 561–604 · score 0.72 · Seurat cluster, FindClusters, FindNeighbors, RunUMAP, harmony, cells
  10. [10] § Methods › Cell type classification › UMAP and Seurat cluster-based cell classification ↔ Rscripts/Neuronal_Glial_Types_Analysis_60D.Rmd, lines 607–644 · score 0.72 · Seurat cluster, FindClusters, FindNeighbors, RunUMAP, harmony, cells
  11. [11] § Methods › Machine learning image models ↔ Image_models/model/data.py, lines 12–34 · score 0.69 · vertical flips, random horizontal, PyTorch, split, batch, training
  12. [12] § Methods › Drug screen on iPSC-derived neural cultures ↔ feature models/sanfilipo_cls_machine_learning.ipynb, lines 924–989 · score 0.68 · Combined antioxidant, Naltrexone, Piracetam, Cannabidiol, Lithium, Carbamazepine
  13. [13] § Methods › Machine learning feature models ↔ feature models/utils/hyperopt.py, lines 26–62 · score 0.66 · cross validation, XGBoost, stratified, Hyperopt, split, fold
  14. [14] § Methods › Cell type classification › Calculation of cell type scores ↔ Rscripts/Neuronal_Glial_Types_Analysis_60D.Rmd, lines 168–228 · score 0.64 · Immature Neurons, AddModuleScore, Pericytes, Nowakowski, OPCs, neuronal
  15. [15] § Results › Machine learning and phenotypic analysis in high-throughput drug screen highlight drugs effective across multiple assays ↔ feature models/sanfilipo_cls_machine_learning.ipynb, lines 924–989 · score 0.63 · machine learning, huperzine, quercetin, resveratrol, dorsomorphin, desferrioxamine
  16. [16] § Results › Machine learning predicts neuroprotective agents effective in high-throughput drug screen ↔ pytorch_grad_cam/seg_eigen_cam.py, lines 12–62 · score 0.59 · Gradient weighted, Grad CAM, neural network, convolutional, exPlanations, Activation
  17. [17] § Methods › Machine learning image models ↔ Image_models/model/optimize.py, lines 15–73 · score 0.59 · Optuna, ADAM, PyTorch, dropout, loss, accuracy
  18. [18] § Results › Machine learning predicts neuroprotective agents effective in high-throughput drug screen ↔ feature models/sanfilipo_cls_machine_learning.ipynb, lines 382–485 · score 0.57 · log odds, Machine learning, feature models, Boxplots, median, boxes
  19. [19] § Methods › Cell type classification › Final cell type assignments ↔ Rscripts/Neuronal_Glial_Types_Analysis_30D.Rmd, lines 606–672 · score 0.55 · cluster_type1, scMayoMap, Seurat cluster, Excitatory, Inhibitory, neurons
  20. [20] § Results › Machine learning predicts neuroprotective agents effective in high-throughput drug screen ↔ feature models/utils/hyperopt.py, lines 26–62 · score 0.55 · cross validation, XGBoost, trees, split, fold, accuracy
  21. [21] § Methods › Cell type classification › Final cell type assignments ↔ Rscripts/Neuronal_Glial_Types_Analysis_60D.Rmd, lines 347–435 · score 0.54 · cluster_type1, scMayoMap, Seurat cluster, Excitatory, Inhibitory, cell
  22. [22] § Methods › Classification of reactive astrocytes ↔ Rscripts/Neuronal_Glial_Types_Analysis_30D.Rmd, lines 731–771 · score 0.51 · Reactive astrocyte cell, AddModuleScore, subset, scores, gene
  23. [23] § Methods › Classification of reactive astrocytes ↔ Rscripts/Neuronal_Glial_Types_Analysis_60D.Rmd, lines 757–798 · score 0.51 · Reactive astrocyte cell, AddModuleScore, subset, scores, gene

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

R Markdown · 940 lines · 56 KB · no license · 5 matches

  1. ---
  2. title: "RNAseq analysis"
  3. author: "Inushi De Silva"
  4. date: '2022-07-12'
  5. output: github_document
  6. editor_options:
  7. chunk_output_type: console
  8. ---
  9. ```{r setup, include=FALSE}
  10. knitr::opts_chunk$set(echo = TRUE)
  11. options("yaml.eval.expr" = TRUE)
  12. ```
  14. 1. ***First we will load the required libraries and set workingdirectories, datadirectories and plotdirectories.***
  15. ```{r Load Libraries and Set Directory Paths}
  16. script_path <- rstudioapi::getActiveDocumentContext()$path
  17. script_dir <- dirname(script_path)
  18. workingdirectory <- dirname(script_dir)
  19. setwd(workingdirectory)
  20. datadirectory <- paste0(workingdirectory,"/data")
  21. resultsdirectory <- paste0(workingdirectory, "/results")
  22. RData_directory <- paste0(datadirectory,"/RData")
  23. ##dir.create(resultsdirectory)
  25. ## Load required libraries
  26. if(!require(dplyr)){
  27. install.packages("dplyr")
  28. library(dplyr)
  29. }
  30. if(!require(Seurat)){
  31. install.packages("Seurat")
  32. library(Seurat)
  33. }
  34. if(!require(ggplot2)){
  35. install.packages("ggplot2")
  36. library(ggplot2)
  37. }
  38. if(!require(sctransform)){
  39. install.packages("sctransform")
  40. library(sctransform)
  41. }
  42. if(!require(RColorBrewer)){
  43. install.packages("RColorBrewer")
  44. library(RColorBrewer)
  45. }
  46. if(!require(ggpointdensity)){
  47. install.packages("ggpointdensity")
  48. library(ggpointdensity)
  49. }
  50. if(!require(devtools)){
  51. install.packages("devtools")
  52. library(devtools)
  53. }
  54. if(!require(BiocManager)){
  55. install.packages("BiocManager")
  56. library(BiocManager)
  57. }
  59. source(file = paste0(script_dir,"/functions_for_sanfil.R"))
  60. ```
  62. 2. ***Load the required RData files that have previously been QC'd and submitted***
  63. ```{r Load RData files that have previously been QC'd and submitted}
  64. ## Load required datasets
  65. set.seed(1)
  66. file_names <- list.files(paste0(RData_directory,"/60D/Sanfil_submitted_versions"), pattern = "*RData",full.names = T)
  67. lapply(file_names[c(5,7)], load, .GlobalEnv)
  68. rm(file_names)
  69. ```
  71. 3. ***We now load in genes from the Nowakowski and Liddelow datasets***
  72. ```{r Load in neuronal gene lists and create colour scheme}
  73. # library(stringi)
  74. # neuronal_genes <- read.delim(file = paste0(datadirectory, "/reference_sheets/neuronal_genelist_new.txt"))
  75. #
  76. # neuronal_genes_ls <- list()
  77. # for(i in 1:ncol(neuronal_genes)){
  78. # neuronal_genes_ls[[i]] <- neuronal_genes[,i]
  79. # names(neuronal_genes_ls)[[i]] <- colnames(neuronal_genes[i])
  80. # }
  81. # neuronal_genes_ls <- rev(neuronal_genes_ls)
  82. #
  83. # neuronal_genes_ls <- lapply(neuronal_genes_ls, function(x){
  84. # stri_remove_empty(x, na_empty = FALSE)
  85. # })
  86. # names(neuronal_genes_ls)[which(names(neuronal_genes_ls) == "Astocyte")] <- "Astrocyte"
  88. ##load updated neuronal gene lists
  89. ##The gene names have been updated below.
  90. load(file = paste0(datadirectory, "/reference_sheets/updated_neuronal_genes_ls.RData"))
  91. neuronal_genes_ls <- updated_neuronal_genes_ls
  93. ## We also create a colour palette for each of the cell types
  94. library(colorspace)
  95. colour_palette <- data.frame(matrix(nrow = 44, ncol = 2))
  96. colour_palette$X1 <- c(names(neuronal_genes_ls)[1:43], "Other")
  97. colour_palette$X2[grepl("Astrocyte", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("salmon", "red"))(3)
  98. colour_palette$X2[grepl("eN", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("limegreen", "palegreen"))(18)
  99. colour_palette$X2[grepl("RG", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("chocolate", "tan"))(11)
  100. colour_palette$X2[grepl("MGE_newborn", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("seagreen", "seagreen1"))(5)
  101. colour_palette$X2[grepl("MGE_Ctx", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("olivedrab", "olivedrab3"))(2)
  102. colour_palette$X2[grepl("CGE_LGE", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("springgreen", "springgreen4"))(2)
  103. colour_palette$X2[grepl("OPC", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("royalblue", "royalblue"))(1)
  104. colour_palette$X2[grepl("Striatal", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("greenyellow", "green"))(1)
  105. colour_palette$X2[grepl("Microglia", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("aquamarine3", "aquamarine3"))(1)
  106. colour_palette$X2[grepl("Mural", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("khaki", "khaki"))(1)
  107. colour_palette$X2[grepl("Other", colour_palette$X1, ignore.case = T)] <- colorRampPalette(c("black", "black"))(1)
  109. ##convert to named vector
  110. colours <- as.vector(colour_palette$X2)
  111. names(colours) <- colour_palette$X1
  112. ```
  114. 4. ***Add genotyping information***
  115. ```{r Add genotyping information - BIAD_ids}
  116. demultiplex_doublet <- read.delim(file = paste0(workingdirectory, "/Cellranger_sheets/60D/Sanfilippo_novaseq/demuxafy/control_60D_combined_vireoSNP/control_60D_combined_results_w_combined_assignments.tsv"))
  117. rownames(demultiplex_doublet) <- demultiplex_doublet$Barcode
  118. cell_names <- rownames([email hidden])
  119. demultiplex_doublet <- demultiplex_doublet[cell_names,]
  120. Control_filtered_ambientRNA_60D <- AddMetaData(Control_filtered_ambientRNA_60D, demultiplex_doublet[,10:11], c("doublets_new","BIAD_id"))
  122. demultiplex_doublet <- read.delim(file = paste0(workingdirectory, "/Cellranger_sheets/60D/Sanfilippo_novaseq/demuxafy/MPSIIIA_60D_combined_vireoSNP/MPSIIIA_combined_results_vireoSNP_w_combined_assignments.tsv"))
  123. rownames(demultiplex_doublet) <- demultiplex_doublet$Barcode
  124. cell_names <- rownames([email hidden])
  125. demultiplex_doublet <- demultiplex_doublet[cell_names,]
  126. MPSIIIA_filtered_ambientRNA_60D <- AddMetaData(MPSIIIA_filtered_ambientRNA_60D, demultiplex_doublet[,10:11], c("doublets_new","BIAD_id"))
  128. DimPlot(Control_filtered_ambientRNA_60D, group.by = "BIAD_id", split.by = "BIAD_id")
  129. table(Control_filtered_ambientRNA_60D$BIAD_id)
  131. save(Control_filtered_ambientRNA_60D, file = paste0(datadirectory, "/RData/60D/revisions/Control_filtered_NOambientRNA_60D.RData"))
  132. save(MPSIIIA_filtered_ambientRNA_60D, file = paste0(datadirectory, "/RData/60D/revisions/MPSIIIA_filtered_NOambientRNA_60D.RData"))
  133. ```
  135. 5. ***Load in processed Lister datasets and markers***
  136. ```{r Load Lister files and markers}
  137. ## Load required datasets
  138. set.seed(1)
  139. file_names <- list.files(paste0(datadirectory,"/Lister_data"), pattern = "*RData",full.names = T)
  140. lapply(file_names[1:2], load, .GlobalEnv)
  141. rm(file_names)
  142. ```
  144. ```{r Prep Lister markers}
  145. # we use genes common between AB and Lister to RunPCA
  146. unique(Lister_markers$cluster) # Astro, OPC, Oligo, Micro,Vas, Inhib: VIP, ID2, PV, SST, PV_SCUBE3, LAMP5_NOS1. Excit: L5-6 THEMIS, L5-6 TLE4, L2-3 CUX2, L4_RORB.
  148. Lister_markers <- Lister_markers[!Lister_markers$cluster %in% "Poor-Quality",]
  149. Lister_markers$cluster <- as.character(Lister_markers$cluster)
  150. Lister_markers$cluster[grep("L4_RORB|L5-6_THEMIS|L2-3_CUX2|L5-6_TLE4",Lister_markers$cluster)] <- "Glutamatergic"
  151. Lister_markers$cluster[grep("LAMP5_NOS1|PV|SST|ID2|VIP",Lister_markers$cluster)] <- "GABAergic"
  152. Lister_markers$cluster[grep("Glut|GABA",Lister_markers$cluster, invert = T)] <- "Non-neuronal"
  153. Lister_markers <- Lister_markers[!Lister_markers$gene[grep("Glut", Lister_markers$cluster)] %in% Lister_markers$gene[grep("GABA", Lister_markers$cluster)],]
  154. Lister_markers <- Lister_markers[!Lister_markers$gene[grep("Glut|GABA", Lister_markers$cluster)] %in% Lister_markers$gene[grep("GABA|Glut", Lister_markers$cluster, invert = T)],]
  156. features_to_use <- process_de_genes(Lister_markers,500,log2FC = 1,p_val_adj = 0.05)
  157. ```
  159. XX. ***Load the required RData files that are being revised***
  160. ```{r Load RData files that have previously been QC'd}
  161. ## Load required datasets
  162. set.seed(1)
  163. file_names <- list.files(paste0(RData_directory,"/60D/revisions"), pattern = "*RData",full.names = T)
  164. lapply(file_names[c(4,5,6,7,8,9)], load, .GlobalEnv)
  165. rm(file_names)
  166. ```
  168. ## ----- ***Analyse the Healthy 60D neurons***
  169. 1. ***We use the Nowakowski and Liddelow dataset to determine cell types in the Control***
  170. ```{r Using the Nowakowski and Liddelow dataset to determine cell types - Control}
  171. Control_filtered_ambientRNA_60D <- AddModuleScore(
  172. object = Control_filtered_ambientRNA_60D,
  173. features = neuronal_genes_ls,
  174. name = names(neuronal_genes_ls),
  175. nbin = 30,
  176. ctrl = 100,
  177. seed = 1,
  178. search = TRUE
  179. )
  180. names([email hidden])[16:24] <- gsub('.{1}$', '', names([email hidden])[16:24])
  181. names([email hidden])[25:58] <- gsub('.{2}$', '', names([email hidden])[25:58])
  183. [email hidden]$maxscore_cell_types <- apply([email hidden][,16:58], 1, function(x)x[maxn(1)(x)])
  184. [email hidden]$cell_types <- apply([email hidden][,16:58], 1, function(x) names(x)[maxn(1)(x)])
  185. ## Add threshold to cell_types which have a low score
  186. Control_filtered_ambientRNA_60D$new_cell_types <- ifelse(Control_filtered_ambientRNA_60D$maxscore_cell_types < 0.1, "Other", Control_filtered_ambientRNA_60D$cell_types)
  189. [email hidden]$maxscore_celltypes_nowakowski <- apply([email hidden][,18:58], 1, function(x)x[maxn(1)(x)])
  190. [email hidden]$nowakowski_celltypes <- apply([email hidden][,18:58], 1, function(x) names(x)[maxn(1)(x)])
  192. ## Add threshold to cell_types which have a low score
  193. Control_filtered_ambientRNA_60D$nowakowski_other <- ifelse(Control_filtered_ambientRNA_60D$maxscore_celltypes_nowakowski < 0.1, "Other", Control_filtered_ambientRNA_60D$nowakowski_celltypes)
  195. Control_filtered_ambientRNA_60D$broad_cell_types <- "EMPTY"
  197. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("CGE|MGE|iN|Striatal"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Inhibitory"
  199. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("eN_"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Excitatory"
  201. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("IPC_eN|MGE_progenitors|Newborn"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Immature Neurons"
  203. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("astrocyte"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Astrocyte"
  205. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("RG"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "RG"
  207. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("OPC"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "OPC"
  209. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("endothelial"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Vas"
  211. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("OPC"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "OPC"
  213. Control_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("Pericyte","Other"), collapse = "|"), Control_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Other"
  216. ## We create a new metadata column to identify donors from Control and MPSIIIA
  217. Control_filtered_ambientRNA_60D$condition_donor <- paste(Control_filtered_ambientRNA_60D$orig.ident, Control_filtered_ambientRNA_60D$donor_id, sep = "_")
  218. Control_filtered_ambientRNA_60D$condition_broadtypes <- paste(Control_filtered_ambientRNA_60D$orig.ident, Control_filtered_ambientRNA_60D$broad_cell_types, sep = "_")
  220. tbl_Control <- [email hidden] %>%
  221. group_by(BIAD_id) %>%
  222. summarise(n = n()) %>%
  223. mutate(freq = (n / sum(n))*100)
  224. tbl_Control
  225. write.csv(tbl_Control, file = paste0(resultsdirectory, "/sheets/revisions_60D/celltype_proportions_60D/Healthy_BIAD_ID_prop.csv"))
  227. save(Control_filtered_ambientRNA_60D, file = paste0(RData_directory, "/60D/revision/Control_filtered_NOambientRNA_60D.RData"))
  228. ```
  230. 2. ***Use Tirosh dataset to determine the cell cycle***
  231. ```{r Using the Tirosh dataset to determine cell cycle - Control}
  232. cell_cycle <- read.table(file = paste0(datadirectory,"/reference_sheets/cellcycle_tirosh.tsv"), header = 1)
  233. cell_cycle <- list("G1S" = cell_cycle$G1.S, "G2M" = cell_cycle$G2.M)
  234. cell_cycle$G1S <- na.omit(cell_cycle$G1S)
  236. Control_filtered_ambientRNA_60D <- AddModuleScore(
  237. object = Control_filtered_ambientRNA_60D,
  238. features = cell_cycle,
  239. name = names(cell_cycle),
  240. nbin = 30,
  241. ctrl = 100,
  242. seed = 1,
  243. search = TRUE
  244. )
  245. names([email hidden])[72:73] <- c("G1S","G2M")
  246. [email hidden]$maxscore_cell_cycle <- apply([email hidden][,72:73], 1, function(x)x[maxn(1)(x)])
  247. [email hidden]$cell_cycle <- apply([email hidden][,72:73], 1, function(x) names(x)[maxn(1)(x)])
  249. Idents(Control_filtered_ambientRNA_60D) <- "broad_cell_types"
  250. DimPlot(Control_filtered_ambientRNA_60D, split.by = "cell_cycle", pt.size = 2)
  252. control_matrix <- as.matrix(Control_filtered_ambientRNA_60D@assays$SCT@data)
  253. control_matrix <- lapply(cell_cycle, function(x){
  254. data.frame(t(as.matrix(control_matrix[which(rownames(control_matrix) %in% x),])))})
  256. for(i in 1:length(control_matrix)){
  257. control_matrix[[i]]$Genes_expressed <- rowSums(control_matrix[[i]] > 0)
  258. }
  260. Control_filtered_ambientRNA_60D <- AddMetaData(Control_filtered_ambientRNA_60D,((control_matrix$G1S$Genes_expressed/43)*100),"G1S_percent_genes_expressed")
  262. Control_filtered_ambientRNA_60D <- AddMetaData(Control_filtered_ambientRNA_60D,((control_matrix$G2M$Genes_expressed/53)*100),"G2M_percent_genes_expressed")
  264. Control_filtered_ambientRNA_60D$percent_genes_expressed <- ifelse(Control_filtered_ambientRNA_60D$cell_cycle == "G1S", Control_filtered_ambientRNA_60D$G1S_percent_genes_expressed, Control_filtered_ambientRNA_60D$G2M_percent_genes_expressed)
  266. metadata <- [email hidden][,c(64,72:78)]
  267. metadata <- split(metadata, metadata$cell_cycle)
  269. library(scater)
  270. ggplot <- lapply(metadata, function(x){
  271. ggplot(x, aes(maxscore_cell_cycle,
  272. percent_genes_expressed, color=broad_cell_types)) +
  273. ##geom_point(aes(color = x$MaxScore_All), alpha = 0.5, size = 1) +
  274. geom_pointdensity(adjust = 0.2) +
  275. geom_point(alpha = 0.5, size = 1) +
  276. theme_minimal() +
  277. xlab("Cell cycle score") + ylab("Percent genes expressed") +
  278. # scale_color_manual(values = c("#80CDC1", ## LSpp
  279. # "#018571", ## LSppCVH
  280. # "#DFC27D", ## TLpp
  281. # "#A6611A")) +
  282. ##scale_colour_gradientn(colours = colorRampPalette(rev(brewer.pal(7,'RdBu')))(9)) +
  283. ##scale_colour_gradientn(colours = pals::coolwarm(7),limits = c(-0.1,0.1), na.value = "grey89") +
  284. geom_vline(xintercept = c(0,0.1), linetype = "dashed", col = "black", linewidth = 0.5) +
  285. geom_hline(yintercept = 30, linetype = "dashed", col = "black", linewidth = 0.5) +
  286. theme(legend.position = "bottom", axis.text.x = element_text(size = 12, hjust = 0.5, angle = 0, vjust = 0.5), axis.text.y = element_text(size = 12), axis.title.x = element_text(size = 12), axis.title.y = element_text(size = 12), panel.border = element_rect(fill=NA, colour = "black", size=1.5)) +
  287. scale_x_continuous(breaks = seq(-0.1,2, by = 0.5)) +
  288. scale_y_continuous(breaks = seq(0, 100 , by = 25)) +
  289. coord_cartesian(ylim = c(0, 100), xlim = c(-0.1,2))
  290. })
  291. gridExtra::grid.arrange(grobs = ggplot, ncol = 2)
  293. [email hidden]$cell_cycle <- ifelse(Control_filtered_ambientRNA_60D$cell_cycle == "G2M" & Control_filtered_ambientRNA_60D$maxscore_cell_cycle < 0.1 & Control_filtered_ambientRNA_60D$percent_genes_expressed <= 25, "NonCycling", ifelse(Control_filtered_ambientRNA_60D$cell_cycle == "G1S" & Control_filtered_ambientRNA_60D$maxscore_cell_cycle < 0.1 & Control_filtered_ambientRNA_60D$percent_genes_expressed <= 25, "NonCycling", Control_filtered_ambientRNA_60D$cell_cycle))
  295. Idents(Control_filtered_ambientRNA_60D) <- "broad_cell_types"
  296. DimPlot(Control_filtered_ambientRNA_60D, split.by = "cell_cycle", pt.size = 2)
  297. ```
  299. 3. ***Integrate with Lister biopsies***
  300. ```{r Integrate 120-day datasets with AB cortical datasets}
  301. genes_to_keep <- Reduce(intersect, list(rownames(Control_filtered_ambientRNA_60D@assays$RNA$counts), rownames(lister_seurat@assays$RNA$counts)))
  303. Combined_Lister_hiPSC_Control_60D <- merge(lister_seurat, c(Control_filtered_ambientRNA_60D))
  304. DefaultAssay(Combined_Lister_hiPSC_Control_60D) <- "RNA"
  305. Combined_Lister_hiPSC_Control_60D[["SCT"]] <- NULL
  306. Combined_Lister_hiPSC_Control_60D[["RNA"]]$data.2 <- NULL
  308. Combined_Lister_hiPSC_Control_60D <- subset(Combined_Lister_hiPSC_Control_60D, features = genes_to_keep)
  310. Combined_Lister_hiPSC_Control_60D$orig.ident[Combined_Lister_hiPSC_Control_60D$orig.ident == "Control"] <- "hiPSC"
  311. Combined_Lister_hiPSC_Control_60D$orig.ident[Combined_Lister_hiPSC_Control_60D$orig.ident == "SeuratProject"] <- "Lister"
  313. Combined_Lister_hiPSC_Control_60D <- SCTransform(Combined_Lister_hiPSC_Control_60D, assay = "RNA")
  315. pc_features <- rownames(Combined_Lister_hiPSC_Control_60D@assays$SCT$scale.data)[rownames(Combined_Lister_hiPSC_Control_60D@assays$SCT$scale.data) %in% unlist(unname(features_to_use))]
  317. Combined_Lister_hiPSC_Control_60D <- RunPCA(Combined_Lister_hiPSC_Control_60D, features = pc_features)
  318. min.pc <- significant_pcs(Combined_Lister_hiPSC_Control_60D,"pca")
  320. #load(file = paste0(datadirectory,"/RData/Biopsy_RData/multicortex_hiPSC_combined_qUMI.RData"))
  321. gc()
  322. library(harmony)
  323. Combined_Lister_hiPSC_Control_60D_int <- RunHarmony(Combined_Lister_hiPSC_Control_60D, "orig.ident", reduction.use = "pca", reduction.save = "harmony", early_stop =F, max_iter = 20, ncores = 4, plot_convergence = T)
  325. min.pc <- significant_pcs(Combined_Lister_hiPSC_Control_60D_int,"harmony")
  327. Combined_Lister_hiPSC_Control_60D_int <- FindNeighbors(Combined_Lister_hiPSC_Control_60D_int, dims = 1:min.pc, reduction = "harmony")
  328. Combined_Lister_hiPSC_Control_60D_int <- FindClusters(Combined_Lister_hiPSC_Control_60D_int, resolution = 0.1)
  329. Combined_Lister_hiPSC_Control_60D_int <- RunUMAP(Combined_Lister_hiPSC_Control_60D_int, dims = 1:min.pc, reduction = "harmony", reduction.name = "umap.int", return.model = T)
  331. #DimPlot(Combined_Lister_hiPSC_Control_60D_int, reduction = "umap.int", group.by = "orig.ident", pt.size = 1, order = F, label = T) +
  332. DimPlot(Combined_Lister_hiPSC_Control_60D_int, reduction = "umap.int", group.by = "new_cell_types", pt.size = 1, order = F, label = T) +
  333. DimPlot(Combined_Lister_hiPSC_Control_60D_int, reduction = "umap.int", group.by = "major_clust", pt.size = 1, order = F, label = T) #+
  334. #DimPlot(Combined_Lister_hiPSC_Control_60D_int, reduction = "umap.int", group.by = "BIAD_id", pt.size = 1, order = F, label = T) & NoLegend()
  336. tbl_Control <- [email hidden] %>%
  337. group_by(BIAD_id) %>%
  338. summarise(n = n()) %>%
  339. mutate(freq = (n / sum(n))*100)
  340. tbl_Control
  342. job::job({
  343. save(Combined_Lister_hiPSC_Control_60D,Combined_Lister_hiPSC_Control_60D_int, file = paste0(RData_directory,"/60D/revisions/Control_Lister_Integrated_60D.RData"))
  344. save(Control_filtered_ambientRNA_60D, file = paste0(RData_directory,"/60D/revisions/Control_filtered_NOambientRNA_60D.RData"))})
  345. ```
  347. 4. ***scMayoMap cell annotations and validation***
  348. ```{r scMayo Map cell annotations}
  349. # Find conserved markers across seurat clusters in AB and hiPSC genes
  350. # Put the conserved markers through scMayo
  351. Combined_Lister_hiPSC_Control_60D_int_2 <- PrepSCTFindMarkers(Combined_Lister_hiPSC_Control_60D_int, assay = "SCT")
  352. Idents(Combined_Lister_hiPSC_Control_60D_int_2) <- "seurat_clusters"
  354. conserved <- list()
  355. for(i in unique(Combined_Lister_hiPSC_Control_60D_int_2$seurat_clusters)){
  356. conserved[[i]] <- FindConservedMarkers(Combined_Lister_hiPSC_Control_60D_int_2, ident.1 = i, grouping.var = "orig.ident", min.pct = 0.25, logfc.threshold = 0.5, assay = "SCT", only.pos = T)
  357. }
  358. save(conserved, file = paste0(RData_directory,"/60D/revisions/conserved_markers_Lister_Control_60D.RData"))
  360. load(file = paste0(RData_directory,"/60D/revisions/conserved_markers_Lister_Control_60D.RData"))
  361. conserved_markers <- conserved
  362. for(i in seq_along(conserved_markers)){
  363. if (nrow(conserved_markers[[i]]) > 0) { # Only proceed if the dataframe is not empty
  364. conserved_markers[[i]]$cluster <- names(conserved_markers)[i]
  365. conserved_markers[[i]]$gene <- rownames(conserved_markers[[i]])
  366. }
  367. }
  369. for(i in 1:length(conserved_markers)){
  370. conserved_markers[[i]] <- conserved_markers[[i]] [,-c(grep("AB|Lister|max|min", colnames(conserved_markers[[i]] )))]
  371. colnames(conserved_markers[[i]] ) <- gsub("hiPSC_","", colnames(conserved_markers[[i]] ))
  372. }
  373. conserved_markers <- do.call("rbind", conserved_markers)
  375. library(scMayoMap)
  376. brain_db <- scMayoMapDatabase[grep("brain", scMayoMapDatabase$tissue), c(1:2,grep("brain", colnames(scMayoMapDatabase)))]
  377. brain_db <- brain_db[,colnames(brain_db)[c(1:2,5,7,10:11,16:17,19,21:22,26,29,31:33,39:40,45:46)]]
  378. rownames(brain_db) <- brain_db$gene
  379. brain_db <- as.data.frame(t(brain_db))
  380. brain_db$cell_type <- gsub("brain:","",rownames(brain_db))
  382. brain_db <- split(brain_db, brain)
  383. for(i in names(brain_db)){
  384. rownames(brain_db[[i]]) <- brain_db[[i]]$gene
  385. brain_db[[i]] <- brain_db[[i]][,names(brain_db[[i]]) %in% i, drop = F]
  386. }
  388. scMayo_cluster <- scMayoMap(conserved_markers, padj.cutoff = 0.05, database = brain_db, pct.cutoff = 0.25)
  389. scMayoplot <- scMayoMap.plot(scMayo_cluster) + scale_colour_gradient(low = "grey70", high = "red") + theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))
  390. ggsave(scMayoplot, filename = paste0(resultsdirectory,"/plots/revisions_60D/[DOTPLOT] Healthy_scMayo_predictions.pdf", sep = ""), height = 6.5, width = 12, dpi = 300, useDingbats = FALSE, limitsize = F)
  392. scMayoMap.plot(scMayo_cluster) + DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "seurat_clusters", pt.size = 0.5, label = T) + DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "major_clust", pt.size = 0.5, label = T) & NoLegend()
  394. scMayo.norm <- scMayo_cluster$annotation.norm[,grep("brain",colnames(scMayo_cluster$annotation.norm))]
  395. scMayo.norm <- scMayo.norm[,colSums(scMayo.norm) > 0]
  396. scMayo.norm$max <- apply(scMayo.norm, 1, function(x)names(x)[maxn(1)(x)]) # highest value
  397. scMayo.norm$max_2 <- apply(scMayo.norm[,1:17], 1, function(x)names(x)[maxn(2)(x)]) #second highest value
  398. scMayo.norm$max <- gsub("brain:","",scMayo.norm$max)
  399. scMayo.norm$max_2 <- gsub("brain:","",scMayo.norm$max_2)
  400. scMayo.norm$integrated_clusters <- rownames(scMayo.norm)
  402. Combined_Lister_hiPSC_Control_60D_int$cluster_type1 <- NA
  403. Combined_Lister_hiPSC_Control_60D_int$cluster_type2 <- NA
  405. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- scMayo.norm$max[match(Combined_Lister_hiPSC_Control_60D_int$seurat_clusters[Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"], scMayo.norm$integrated_clusters)]
  407. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- scMayo.norm$max_2[match(Combined_Lister_hiPSC_Control_60D_int$seurat_clusters[Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"], scMayo.norm$integrated_clusters)]
  409. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("Glut",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "Excitatory"
  410. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("Glut",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "Excitatory"
  412. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("GABA|Interneuron|basket|Basket|Chandelier|Martinotti",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "Inhibitory"
  413. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("GABA|Interneuron|basket|Basket|Chandelier|Martinotti",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "Inhibitory"
  415. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("precursor",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "OPC"
  416. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("precursor",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "OPC"
  418. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("Oligodendrocyte",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "Oligo"
  419. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("Oligodendrocyte",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "Oligo"
  421. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("Radial",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "RG"
  422. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("Radial",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "RG"
  424. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("progenitor|stem|Neuroblast",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "RG"
  425. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("progenitor|stem|Neuroblast",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "RG"
  427. Combined_Lister_hiPSC_Control_60D_int$cluster_type1[grep("Micro",Combined_Lister_hiPSC_Control_60D_int$cluster_type1)] <- "Micro"
  429. Combined_Lister_hiPSC_Control_60D_int$cluster_type2[grep("Endo",Combined_Lister_hiPSC_Control_60D_int$cluster_type2)] <- "Vas"
  431. DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "cluster_type1", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "cluster_type2", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "major_clust", label = T, order = T) #& NoLegend()
  433. job::job({
  434. save(Combined_Lister_hiPSC_Control_60D,Combined_Lister_hiPSC_Control_60D_int, file = paste0(RData_directory,"/60D/revisions/Control_Lister_Integrated_60D.RData"))})
  435. ```
  437. ```{r Match assignments between Nowakowski, Lister and scMayo}
  438. DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "major_clust", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "cluster_type1", label = T, order = T) +DimPlot(Combined_Lister_hiPSC_Control_60D_int, group.by = "seurat_clusters", label = T, order = T)#& NoLegend()
  440. tbl_Control <- [email hidden] %>%
  441. group_by(seurat_clusters, major_clust) %>%
  442. summarise(n = n()) %>%
  443. mutate(freq = (n / sum(n))*100)
  444. tbl_Control
  446. Combined_Lister_hiPSC_Control_60D_int$cell_type[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(4,6,7,8) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Inhibitory"
  447. Combined_Lister_hiPSC_Control_60D_int$cell_type[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(2,10,14,12) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Excitatory"
  448. Combined_Lister_hiPSC_Control_60D_int$cell_type[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(0,1,3,5,9,11,13) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Non-neuronal"
  450. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(4,6,7,8,12) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Inhibitory"
  451. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(2,10,14,12) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Excitatory"
  452. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(1,9) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Astrocyte"
  453. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(5) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Micro"
  454. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(0) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Oligo"
  455. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(3) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "OPC"
  456. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$seurat_clusters %in% c(13,11) & Combined_Lister_hiPSC_Control_60D_int$orig.ident == "hiPSC"] <- "Vas"
  458. Combined_Lister_hiPSC_Control_60D_int$major_clust[Combined_Lister_hiPSC_Control_60D_int$major_clust == "Astro"] <- "Astrocyte"
  460. metadata <- [email hidden][[email hidden]$orig.ident == "hiPSC",c("major_clust", "cluster_type1","cluster_type2", "broad_cell_types")]
  462. matched_cells <- assign_final_cell_type(metadata,majority_cols = c("major_clust", "cluster_type1","cluster_type2", "broad_cell_types"),major_clust = "major_clust", tiebreaker_col = "major_clust", ignore_values = c("Vas","Oligo","Micro"), mix_col = "broad_cell_types")
  464. matched_cells$majority_match[is.na(matched_cells$majority_match)] <- ifelse(matched_cells$major_clust[is.na(matched_cells$majority_match)] == "Astrocyte" & matched_cells$cluster_type2[is.na(matched_cells$majority_match)] == "RG" & matched_cells$broad_cell_types[is.na(matched_cells$majority_match)] %in% c("RG","OPC","Other"), "RG","Immature Neurons")
  466. Control_filtered_ambientRNA_60D <- AddMetaData(Control_filtered_ambientRNA_60D, matched_cells[,6], c("final_cell_type"))
  467. DimPlot(Control_filtered_ambientRNA_60D, group.by = "final_cell_type",pt.size = 1)
  469. table(Control_filtered_ambientRNA_60D$final_cell_type)
  471. tbl_Control <- [email hidden]%>%
  472. group_by(BIAD_id,final_cell_type,cell_cycle) %>%
  473. summarise(n = n()) %>%
  474. mutate(freq = (n / sum(n))*100)
  475. tbl_Control
  476. write.csv(tbl_Control, file = paste0(resultsdirectory, "/sheets/revisions_60D/celltype_proportions_60D/Healthy_BIAD_ID_cell_type_cycling_prop.csv"))
  478. save(Control_filtered_ambientRNA_60D, file = paste0(RData_directory,"/60D/revisions/Control_filtered_NOambientRNA_60D.RData"))
  479. ```
  481. ## ----- ***Analyse the MPSIIIA 60D neurons***
  483. 1. ***We use the Nowakowski and Liddelow dataset to determine cell types in Sanfilippo***
  484. ```{r Using the Nowakowski and Liddelow dataset to determine cell types - Sanfilippo}
  485. MPSIIIA_filtered_ambientRNA_60D <- AddModuleScore(
  486. object = MPSIIIA_filtered_ambientRNA_60D,
  487. features = neuronal_genes_ls,
  488. name = names(neuronal_genes_ls),
  489. nbin = 30,
  490. ctrl = 100,
  491. seed = 1,
  492. search = TRUE
  493. )
  494. names([email hidden])[16:24] <- gsub('.{1}$', '', names([email hidden])[16:24])
  495. names([email hidden])[25:58] <- gsub('.{2}$', '', names([email hidden])[25:58])
  497. [email hidden]$maxscore_cell_types <- apply([email hidden][,16:58], 1, function(x)x[maxn(1)(x)])
  498. [email hidden]$cell_types <- apply([email hidden][,16:58], 1, function(x) names(x)[maxn(1)(x)])
  499. ## Add threshold to cell_types which have a low score
  500. MPSIIIA_filtered_ambientRNA_60D$new_cell_types <- ifelse(MPSIIIA_filtered_ambientRNA_60D$maxscore_cell_types < 0.1, "Other", MPSIIIA_filtered_ambientRNA_60D$cell_types)
  502. [email hidden]$maxscore_celltypes_nowakowski <- apply([email hidden][,18:58], 1, function(x)x[maxn(1)(x)])
  503. [email hidden]$nowakowski_celltypes <- apply([email hidden][,18:58], 1, function(x) names(x)[maxn(1)(x)])
  504. ## Add threshold to cell_types which have a low score
  505. MPSIIIA_filtered_ambientRNA_60D$nowakowski_other <- ifelse(MPSIIIA_filtered_ambientRNA_60D$maxscore_celltypes_nowakowski < 0.1, "Other", MPSIIIA_filtered_ambientRNA_60D$nowakowski_celltypes)
  507. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types <- "EMPTY"
  509. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("CGE|MGE|iN|Striatal"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Inhibitory"
  511. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("eN_"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Excitatory"
  513. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("IPC_eN|MGE_progenitors|Newborn"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Immature Neurons"
  515. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("astrocyte"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Astrocyte"
  517. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("RG"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "RG"
  519. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("OPC"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "OPC"
  521. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("endothelial"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Vas"
  523. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("OPC"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "OPC"
  525. MPSIIIA_filtered_ambientRNA_60D$broad_cell_types[grep(paste(c("Pericyte","Other"), collapse = "|"), MPSIIIA_filtered_ambientRNA_60D$new_cell_types, ignore.case = T)] <- "Other"
  527. ## We create a new metadata column to identify donors from Control and MPSIIIA
  528. MPSIIIA_filtered_ambientRNA_60D$condition_donor <- paste(MPSIIIA_filtered_ambientRNA_60D$orig.ident, MPSIIIA_filtered_ambientRNA_60D$donor_id, sep = "_")
  529. MPSIIIA_filtered_ambientRNA_60D$condition_broadtypes <- paste(MPSIIIA_filtered_ambientRNA_60D$orig.ident, MPSIIIA_filtered_ambientRNA_60D$broad_cell_types, sep = "_")
  531. tbl_MPSIIIA <- [email hidden] %>%
  532. group_by(donor_id,BIAD_id) %>%
  533. summarise(n = n()) %>%
  534. mutate(freq = (n / sum(n))*100)
  535. tbl_MPSIIIA
  536. write.csv(tbl_MPSIIIA, file = paste0(resultsdirectory, "/sheets/revisions_60D/celltype_proportions_60D/MPSIIIA_BIAD_ID_w_old_donor_id_prop.csv"))
  538. save(MPSIIIA_filtered_ambientRNA_60D, file = paste0(datadirectory, "/RData/60D/revisions/MPSIIIA_filtered_NOambientRNA_60D.RData"))
  539. ```
  541. 2. ***Use Tirosh dataset to determine the cell cycle***
  542. ```{r Using the Tirosh dataset to determine cell cycle - MPSIIIA}
  543. cell_cycle <- read.table(file = paste0(datadirectory,"/reference_sheets/cellcycle_tirosh.tsv"), header = 1)
  544. cell_cycle <- list("G1S" = cell_cycle$G1.S, "G2M" = cell_cycle$G2.M)
  545. cell_cycle$G1S <- na.omit(cell_cycle$G1S)
  547. MPSIIIA_filtered_ambientRNA_60D <- AddModuleScore(
  548. object = MPSIIIA_filtered_ambientRNA_60D,
  549. features = cell_cycle,
  550. name = names(cell_cycle),
  551. nbin = 30,
  552. ctrl = 100,
  553. seed = 1,
  554. search = TRUE
  555. )
  556. names([email hidden])[72:73] <- c("G1S","G2M")
  557. [email hidden]$maxscore_cell_cycle <- apply([email hidden][,72:73], 1, function(x)x[maxn(1)(x)])
  558. [email hidden]$cell_cycle <- apply([email hidden][,72:73], 1, function(x) names(x)[maxn(1)(x)])
  560. MPSIIIA_matrix <- as.matrix(MPSIIIA_filtered_ambientRNA_60D@assays$SCT@data)
  561. MPSIIIA_matrix <- lapply(cell_cycle, function(x){
  562. data.frame(t(as.matrix(MPSIIIA_matrix[which(rownames(MPSIIIA_matrix) %in% x),])))})
  564. for(i in 1:length(MPSIIIA_matrix)){
  565. MPSIIIA_matrix[[i]]$Genes_expressed <- rowSums(MPSIIIA_matrix[[i]] > 0)
  566. }
  568. MPSIIIA_filtered_ambientRNA_60D <- AddMetaData(MPSIIIA_filtered_ambientRNA_60D,((MPSIIIA_matrix$G1S$Genes_expressed/43)*100),"G1S_percent_genes_expressed")
  570. MPSIIIA_filtered_ambientRNA_60D <- AddMetaData(MPSIIIA_filtered_ambientRNA_60D,((MPSIIIA_matrix$G2M$Genes_expressed/53)*100),"G2M_percent_genes_expressed")
  572. MPSIIIA_filtered_ambientRNA_60D$percent_genes_expressed <- ifelse(MPSIIIA_filtered_ambientRNA_60D$cell_cycle == "G1S", MPSIIIA_filtered_ambientRNA_60D$G1S_percent_genes_expressed, MPSIIIA_filtered_ambientRNA_60D$G2M_percent_genes_expressed)
  574. metadata <- [email hidden][,c(64,72:78)]
  575. metadata <- split(metadata, metadata$cell_cycle)
  577. library(scater)
  578. ggplot <- lapply(metadata, function(x){
  579. ggplot(x, aes(maxscore_cell_cycle,
  580. percent_genes_expressed, color=broad_cell_types)) +
  581. ##geom_point(aes(color = x$MaxScore_All), alpha = 0.5, size = 1) +
  582. geom_pointdensity(adjust = 0.2) +
  583. geom_point(alpha = 0.5, size = 1) +
  584. theme_minimal() +
  585. xlab("Cell cycle score") + ylab("Percent genes expressed") +
  586. # scale_color_manual(values = c("#80CDC1", ## LSpp
  587. # "#018571", ## LSppCVH
  588. # "#DFC27D", ## TLpp
  589. # "#A6611A")) +
  590. ##scale_colour_gradientn(colours = colorRampPalette(rev(brewer.pal(7,'RdBu')))(9)) +
  591. ##scale_colour_gradientn(colours = pals::coolwarm(7),limits = c(-0.1,0.1), na.value = "grey89") +
  592. geom_vline(xintercept = c(0,0.1), linetype = "dashed", col = "black", linewidth = 0.5) +
  593. geom_hline(yintercept = 30, linetype = "dashed", col = "black", linewidth = 0.5) +
  594. theme(legend.position = "bottom", axis.text.x = element_text(size = 12, hjust = 0.5, angle = 0, vjust = 0.5), axis.text.y = element_text(size = 12), axis.title.x = element_text(size = 12), axis.title.y = element_text(size = 12), panel.border = element_rect(fill=NA, colour = "black", size=1.5)) +
  595. scale_x_continuous(breaks = seq(-0.1,2, by = 0.5)) +
  596. scale_y_continuous(breaks = seq(0, 100 , by = 25)) +
  597. coord_cartesian(ylim = c(0, 100), xlim = c(-0.1,2))
  598. })
  599. gridExtra::grid.arrange(grobs = ggplot, ncol = 2)
  601. [email hidden]$cell_cycle <- ifelse(MPSIIIA_filtered_ambientRNA_60D$cell_cycle == "G2M" & MPSIIIA_filtered_ambientRNA_60D$maxscore_cell_cycle < 0.1 & MPSIIIA_filtered_ambientRNA_60D$percent_genes_expressed <= 25, "NonCycling", ifelse(MPSIIIA_filtered_ambientRNA_60D$cell_cycle == "G1S" & MPSIIIA_filtered_ambientRNA_60D$maxscore_cell_cycle < 0.1 & MPSIIIA_filtered_ambientRNA_60D$percent_genes_expressed <= 25, "NonCycling", MPSIIIA_filtered_ambientRNA_60D$cell_cycle))
  603. Idents(MPSIIIA_filtered_ambientRNA_60D) <- "broad_cell_types"
  604. DimPlot(MPSIIIA_filtered_ambientRNA_60D, split.by = "cell_cycle", pt.size = 2)
  605. ```
  607. 3. ***Integrate with Lister biopsies***
  608. ```{r Integrate 120-day datasets with AB cortical datasets}
  609. genes_to_keep <- Reduce(intersect, list(rownames(MPSIIIA_filtered_ambientRNA_60D@assays$RNA$counts), rownames(lister_seurat@assays$RNA$counts)))
  611. Combined_Lister_hiPSC_MPSIIIA_60D <- merge(lister_seurat, c(MPSIIIA_filtered_ambientRNA_60D))
  612. DefaultAssay(Combined_Lister_hiPSC_MPSIIIA_60D) <- "RNA"
  613. Combined_Lister_hiPSC_MPSIIIA_60D[["SCT"]] <- NULL
  614. Combined_Lister_hiPSC_MPSIIIA_60D[["RNA"]]$data.Sanfilippo_60D <- NULL
  616. Combined_Lister_hiPSC_MPSIIIA_60D <- subset(Combined_Lister_hiPSC_MPSIIIA_60D, features = genes_to_keep)
  618. Combined_Lister_hiPSC_MPSIIIA_60D$orig.ident[Combined_Lister_hiPSC_MPSIIIA_60D$orig.ident == "Sanfilippo"] <- "hiPSC"
  619. Combined_Lister_hiPSC_MPSIIIA_60D$orig.ident[Combined_Lister_hiPSC_MPSIIIA_60D$orig.ident == "SeuratProject"] <- "Lister"
  621. Combined_Lister_hiPSC_MPSIIIA_60D <- SCTransform(Combined_Lister_hiPSC_MPSIIIA_60D, assay = "RNA")
  623. pc_features <- rownames(Combined_Lister_hiPSC_MPSIIIA_60D@assays$SCT$scale.data)[rownames(Combined_Lister_hiPSC_MPSIIIA_60D@assays$SCT$scale.data) %in% unlist(unname(features_to_use))]
  625. Combined_Lister_hiPSC_MPSIIIA_60D <- RunPCA(Combined_Lister_hiPSC_MPSIIIA_60D, features = pc_features)
  626. min.pc <- significant_pcs(Combined_Lister_hiPSC_MPSIIIA_60D,"pca")
  628. #load(file = paste0(datadirectory,"/RData/Biopsy_RData/multicortex_hiPSC_combined_qUMI.RData"))
  629. gc()
  630. library(harmony)
  631. Combined_Lister_hiPSC_MPSIIIA_60D_int <- RunHarmony(Combined_Lister_hiPSC_MPSIIIA_60D, "orig.ident", reduction.use = "pca", reduction.save = "harmony", early_stop =F, max_iter = 20, ncores = 4, plot_convergence = T)
  633. min.pc <- significant_pcs(Combined_Lister_hiPSC_MPSIIIA_60D_int,"harmony")
  635. Combined_Lister_hiPSC_MPSIIIA_60D_int <- FindNeighbors(Combined_Lister_hiPSC_MPSIIIA_60D_int, dims = 1:min.pc, reduction = "harmony")
  636. Combined_Lister_hiPSC_MPSIIIA_60D_int <- FindClusters(Combined_Lister_hiPSC_MPSIIIA_60D_int, resolution = 0.1)
  637. Combined_Lister_hiPSC_MPSIIIA_60D_int <- RunUMAP(Combined_Lister_hiPSC_MPSIIIA_60D_int, dims = 1:min.pc, reduction = "harmony", reduction.name = "umap.int", return.model = T)
  639. DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, reduction = "umap.int", group.by = "orig.ident", pt.size = 1, order = F, label = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, reduction = "umap.int", group.by = "major_clust", pt.size = 1, order = F, label = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, reduction = "umap.int", group.by = "broad_cell_types", pt.size = 1, order = F, label = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, reduction = "umap.int", group.by = "seurat_clusters", pt.size = 1, order = F, label = T) & NoLegend()
  641. job::job({
  642. save(Combined_Lister_hiPSC_MPSIIIA_60D,Combined_Lister_hiPSC_MPSIIIA_60D_int, file = paste0(RData_directory,"/60D/revisions/MPSIIIA_Lister_Integrated_60D.RData"))
  643. save(MPSIIIA_filtered_ambientRNA_60D, file = paste0(RData_directory,"/60D/revisions/MPSIIIA_filtered_NOambientRNA_60D.RData"))})
  644. ```
  646. 4. ***scMayoMap cell annotations and validation***
  647. ```{r scMayo Map cell annotations}
  648. # Find conserved markers across seurat clusters in AB and hiPSC genes
  649. # Put the conserved markers through scMayo
  650. Combined_Lister_hiPSC_MPSIIIA_60D_int_2 <- PrepSCTFindMarkers(Combined_Lister_hiPSC_MPSIIIA_60D_int, assay = "SCT")
  651. Idents(Combined_Lister_hiPSC_MPSIIIA_60D_int_2) <- "seurat_clusters"
  653. conserved <- list()
  654. for(i in unique(Combined_Lister_hiPSC_MPSIIIA_60D_int_2$seurat_clusters)){
  655. conserved[[i]] <- FindConservedMarkers(Combined_Lister_hiPSC_MPSIIIA_60D_int_2, ident.1 = i, grouping.var = "orig.ident", min.pct = 0.25, logfc.threshold = 0.5, assay = "SCT", only.pos = T)
  656. }
  657. save(conserved, file = paste0(RData_directory,"/60D/revisions/conserved_markers_Lister_MPSIIIA_60D.RData"))
  659. load(file = paste0(RData_directory,"/60D/revisions/conserved_markers_Lister_MPSIIIA_60D.RData"))
  660. conserved_markers <- conserved
  661. for(i in seq_along(conserved_markers)){
  662. if (nrow(conserved_markers[[i]]) > 0) { # Only proceed if the dataframe is not empty
  663. conserved_markers[[i]]$cluster <- names(conserved_markers)[i]
  664. conserved_markers[[i]]$gene <- rownames(conserved_markers[[i]])
  665. }
  666. }
  668. for(i in 1:length(conserved_markers)){
  669. conserved_markers[[i]] <- conserved_markers[[i]] [,-c(grep("AB|Lister|max|min", colnames(conserved_markers[[i]] )))]
  670. colnames(conserved_markers[[i]] ) <- gsub("hiPSC_","", colnames(conserved_markers[[i]] ))
  671. }
  672. conserved_markers <- do.call("rbind", conserved_markers)
  674. library(scMayoMap)
  675. brain_db <- scMayoMapDatabase[grep("brain", scMayoMapDatabase$tissue), c(1:2,grep("brain", colnames(scMayoMapDatabase)))]
  677. brain_db <- brain_db[,colnames(brain_db)[c(1:2,5,7,10:11,16:17,19,21:22,26,29,31:33,35,39:40,45:46)]]
  678. scMayo_cluster <- scMayoMap(conserved_markers, padj.cutoff = 0.05, database = brain_db, pct.cutoff = 0.25)
  679. scMayoplot <- scMayoMap.plot(scMayo_cluster) + scale_colour_gradient(low = "grey70", high = "red") + theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))
  680. ggsave(scMayoplot, filename = paste0(resultsdirectory,"/plots/revisions_60D/[DOTPLOT] MPSIIIA_scMayo_predictions.pdf", sep = ""), height = 6.5, width = 12, dpi = 300, useDingbats = FALSE, limitsize = F)
  682. scMayoMap.plot(scMayo_cluster) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "seurat_clusters", pt.size = 0.5, label = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "major_clust", pt.size = 0.5, label = T) & NoLegend()
  684. scMayo.norm <- scMayo_cluster$annotation.norm[,grep("brain",colnames(scMayo_cluster$annotation.norm))]
  685. scMayo.norm <- scMayo.norm[,colSums(scMayo.norm) > 0]
  686. scMayo.norm$max <- apply(scMayo.norm, 1, function(x)names(x)[maxn(1)(x)]) # highest value
  687. scMayo.norm$max_2 <- apply(scMayo.norm[,1:18], 1, function(x)names(x)[maxn(2)(x)]) #second highest value
  688. scMayo.norm$max <- gsub("brain:","",scMayo.norm$max)
  689. scMayo.norm$max_2 <- gsub("brain:","",scMayo.norm$max_2)
  690. scMayo.norm$integrated_clusters <- rownames(scMayo.norm)
  692. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1 <- NA
  693. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2 <- NA
  695. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1[Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- scMayo.norm$max[match(Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters[Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"], scMayo.norm$integrated_clusters)]
  697. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2[Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- scMayo.norm$max_2[match(Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters[Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"], scMayo.norm$integrated_clusters)]
  699. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1[grep("Glut|Pyramidal",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1)] <- "Excitatory"
  700. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2[grep("Glut|Pyramidal",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2)] <- "Excitatory"
  702. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1[grep("GABA|Interneuron|Basket|Chandelier|Martinotti",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1)] <- "Inhibitory"
  703. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2[grep("GABA|Interneuron|Basket|Chandelier|Martinotti",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2)] <- "Inhibitory"
  705. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1[grep("precursor",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1)] <- "OPC"
  706. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2[grep("precursor",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2)] <- "OPC"
  708. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1[grep("Oligodendrocyte",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1)] <- "Oligo"
  709. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2[grep("Oligodendrocyte",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2)] <- "Oligo"
  711. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1[grep("Radial|stem",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type1)] <- "RG"
  712. Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2[grep("Radial|stem",Combined_Lister_hiPSC_MPSIIIA_60D_int$cluster_type2)] <- "RG"
  714. DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "cluster_type1", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "cluster_type2", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "major_clust", label = T, order = T) & NoLegend()
  715. ```
  717. ```{r Match assignments between Nowakowski, Lister and scMayo}
  718. DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "major_clust", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "broad_cell_types", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "cluster_type1", label = T, order = T) + DimPlot(Combined_Lister_hiPSC_MPSIIIA_60D_int, group.by = "seurat_clusters", label = T, order = T)& NoLegend()
  720. tbl_MPSIIIA<- [email hidden] %>%
  721. group_by(seurat_clusters, major_clust) %>%
  722. summarise(n = n()) %>%
  723. mutate(freq = (n / sum(n))*100)
  724. tbl_MPSIIIA
  726. Combined_Lister_hiPSC_MPSIIIA_60D_int$cell_type[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(5,7,11) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Inhibitory"
  727. Combined_Lister_hiPSC_MPSIIIA_60D_int$cell_type[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(2,8) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Excitatory"
  728. Combined_Lister_hiPSC_MPSIIIA_60D_int$cell_type[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(0,1,3,4,6,9,10) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Non-neuronal"
  730. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(5,7,11) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Inhibitory"
  731. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(2,8) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Excitatory"
  732. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(1,6) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Astrocyte"
  733. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(4,9) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Micro"
  734. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(0) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Oligo"
  735. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(3) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "OPC"
  736. Combined_Lister_hiPSC_MPSIIIA_60D_int$major_clust[Combined_Lister_hiPSC_MPSIIIA_60D_int$seurat_clusters %in% c(10) & Combined_Lister_hiPSC_MPSIIIA_60D_int$orig.ident == "hiPSC"] <- "Vas"
  738. metadata <- [email hidden][[email hidden]$orig.ident == "hiPSC",c("major_clust", "cluster_type1","cluster_type2", "broad_cell_types")]
  740. matched_cells <- assign_final_cell_type(metadata,majority_cols = c("major_clust", "cluster_type1","cluster_type2", "broad_cell_types"),major_clust = "major_clust", tiebreaker_col = "major_clust", ignore_values = c("Vas","Oligo","Micro"), mix_col = "broad_cell_types")
  742. MPSIIIA_filtered_ambientRNA_60D <- AddMetaData(MPSIIIA_filtered_ambientRNA_60D, matched_cells[,6], c("final_cell_type"))
  743. DimPlot(MPSIIIA_filtered_ambientRNA_60D, group.by = "final_cell_type",pt.size = 1, split.by = "orig.ident")
  745. table(MPSIIIA_filtered_ambientRNA_60D$final_cell_type)
  747. tbl_MPSIIIA<- [email hidden] %>%
  748. group_by(BIAD_id,final_cell_type) %>%
  749. summarise(n = n()) %>%
  750. mutate(freq = (n / sum(n))*100)
  751. tbl_MPSIIIA
  752. write.csv(tbl_MPSIIIA, file = paste0(resultsdirectory, "/sheets/revisions_60D/celltype_proportions_60D/MPSIIIA_BIAD_ID_cell_type_prop.csv"))
  754. save(MPSIIIA_filtered_ambientRNA_60D, file = paste0(RData_directory,"/60D/revisions/MPSIIIA_filtered_NOambientRNA_60D.RData"))
  755. ```
  757. ## ----- **Analyse the proportion of reactive astrocytes**
  758. 1. ***Subset the datasets and calculate astrocyte scores for only the astrocyte cells - Control***
  759. ```{r Subset the datasets and calculate astrocyte scores for only the astrocyte cells - Control}
  760. Idents(Control_filtered_ambientRNA_60D) <- "final_cell_type"
  761. Control_ambientRNA_astro_60D <- subset(Control_filtered_ambientRNA_60D, idents = unique(Control_filtered_ambientRNA_60D$final_cell_type[grep("astrocyte", Control_filtered_ambientRNA_60D$final_cell_type,ignore.case = T)]))
  763. Control_ambientRNA_astro_60D <- SCTransform(Control_ambientRNA_astro_60D)
  764. [email hidden][,76:80] <- NULL
  766. load(file = paste0(datadirectory, "/reference_sheets/updated_neuronal_genes_ls.RData"))
  767. neuronal_genes_ls <- updated_neuronal_genes_ls
  769. Control_ambientRNA_astro_60D <- AddModuleScore(
  770. object = Control_ambientRNA_astro_60D,
  771. features = neuronal_genes_ls[1:2],
  772. name = names(neuronal_genes_ls)[1:2],
  773. nbin = 20, #15 #20
  774. ctrl = 100,
  775. seed = 1,
  776. search = TRUE
  777. )
  778. names([email hidden])[77:78] <- gsub('.{1}$', '', names([email hidden])[77:78])
  780. [email hidden]$maxscore_reactiveastro <- apply([email hidden][,77:78], 1, function(x)x[maxn(1)(x)])
  781. [email hidden]$reactive_astro_type <- apply([email hidden][,77:78], 1, function(x) names(x)[maxn(1)(x)])
  783. Control_ambientRNA_astro_60D$reactive_astro_type_2 <- ifelse(Control_ambientRNA_astro_60D$maxscore_reactiveastro < 0.1 & Control_ambientRNA_astro_60D$final_cell_type != "Astrocyte", "Not Astrocyte", ifelse(Control_ambientRNA_astro_60D$maxscore_reactiveastro < 0.1 & Control_ambientRNA_astro_60D$final_cell_type == "Astrocyte", "Astrocyte", Control_ambientRNA_astro_60D$reactive_astro_type))
  785. [email hidden][,76:80] <- NULL
  787. Control_filtered_ambientRNA_60D <- AddMetaData(Control_filtered_ambientRNA_60D, [email hidden][,77:81], names([email hidden])[77:81])
  789. tbl_Control <- [email hidden] %>%
  790. group_by(BIAD_id, reactive_astro_type_2) %>%
  791. summarise(n = n()) %>%
  792. mutate(freq = (n / sum(n))*100)
  793. tbl_Control
  794. write.csv(tbl_Control, file = paste0(resultsdirectory, "/sheets/revisions_60D/celltype_proportions_60D/Healthy_BIAD_id_reactive_astro_prop.csv"))
  796. save(Control_ambientRNA_astro_60D, file = paste0(datadirectory, "/RData/60D/revisions/Control_NOambientRNA_astro_60D.RData"))
  797. save(Control_filtered_ambientRNA_60D, file = paste0(datadirectory, "/RData/60D/revisions/Control_filtered_NOambientRNA_60D.RData"))
  798. ```
  800. 2. ***Subset the datasets and calculate astrocyte scores for only the astrocyte cells - Sanfilippo***
  801. ```{r Subset the datasets and calculate astrocyte scores for only the astrocyte cells - Sanfilippo}
  802. Idents(MPSIIIA_filtered_ambientRNA_60D) <- "final_cell_type"
  803. MPSIIIA_ambientRNA_astro_60D <- subset(MPSIIIA_filtered_ambientRNA_60D, idents = unique(MPSIIIA_filtered_ambientRNA_60D$final_cell_type[grep("astrocyte", MPSIIIA_filtered_ambientRNA_60D$final_cell_type, ignore.case = T)]))
  805. MPSIIIA_ambientRNA_astro_60D <- SCTransform(MPSIIIA_ambientRNA_astro_60D)
  806. [email hidden][,77:81] <- NULL
  808. MPSIIIA_ambientRNA_astro_60D <- AddModuleScore(
  809. object = MPSIIIA_ambientRNA_astro_60D,
  810. features = neuronal_genes_ls[1:2],
  811. name = names(neuronal_genes_ls)[1:2],
  812. nbin = 20,
  813. ctrl = 100,
  814. seed = 1,
  815. search = TRUE
  816. )
  817. names([email hidden])[77:78] <- gsub('.{1}$', '', names([email hidden])[77:78])
  819. [email hidden]$maxscore_reactiveastro <- apply([email hidden][,77:78], 1, function(x)x[maxn(1)(x)])
  820. [email hidden]$reactive_astro_type <- apply([email hidden][,77:78], 1, function(x) names(x)[maxn(1)(x)])
  822. MPSIIIA_ambientRNA_astro_60D$reactive_astro_type_2 <- ifelse(MPSIIIA_ambientRNA_astro_60D$maxscore_reactiveastro < 0.1 & MPSIIIA_ambientRNA_astro_60D$final_cell_type != "Astrocyte", "Not Astrocyte", ifelse(MPSIIIA_ambientRNA_astro_60D$maxscore_reactiveastro < 0.1 & MPSIIIA_ambientRNA_astro_60D$final_cell_type == "Astrocyte", "Astrocyte", MPSIIIA_ambientRNA_astro_60D$reactive_astro_type))
  825. [email hidden][,77:81] <- NULL
  826. MPSIIIA_filtered_ambientRNA_60D <- AddMetaData(MPSIIIA_filtered_ambientRNA_60D, [email hidden][,77:81], names([email hidden])[77:81])
  828. tbl_MPSIIIA <- [email hidden] %>%
  829. group_by(BIAD_id, reactive_astro_type_2) %>%
  830. summarise(n = n()) %>%
  831. mutate(freq = (n / sum(n))*100)
  832. tbl_MPSIIIA
  833. write.csv(tbl_MPSIIIA, file = paste0(resultsdirectory, "/sheets/revisions_60D/celltype_proportions_60D/MPSIIIA_BIAD_id_reactive_astro_prop.csv"))
  835. save(MPSIIIA_ambientRNA_astro_60D, file = paste0(datadirectory, "/RData/60D/revisions/MPSIIIA_NOambientRNA_astro_60D.RData"))
  837. save(MPSIIIA_filtered_ambientRNA_60D, file = paste0(datadirectory, "/RData/60D/revisions/MPSIIIA_filtered_NOambientRNA_60D.RData"))
  838. ```
  840. ```{r EXTRA - update gene names}
  841. library(openxlsx)
  842. ## Update the gene symbol list - ALL
  843. updated_neuronal_genes_ls <- lapply(neuronal_genes_ls, function(x){
  844. UpdateSymbolList(x,
  845. timeout = 20000)
  846. })
  848. ## reactive astrocytes
  849. updated_neuronal_genes_ls_liddelow <- lapply(neuronal_genes_ls[1:2], function(x){
  850. UpdateSymbolList(x,
  851. timeout = 20000)
  852. })
  854. ## MGE_progenitors
  855. updated_neuronal_genes_ls_nowakowski_MGE_progenitors <- lapply(neuronal_genes_ls[3:6], function(x){
  856. UpdateSymbolList(x,
  857. timeout = 20000)
  858. })
  860. ## MGE RG
  861. updated_neuronal_genes_ls_nowakowski_MGE_RG<- lapply(neuronal_genes_ls[7:8], function(x){
  862. UpdateSymbolList(x,
  863. timeout = 20000)
  864. })
  867. updated_neuronal_genes_ls_nowakowski_DivRG_S<- lapply(neuronal_genes_ls[9], function(x){
  868. UpdateSymbolList(x,
  869. timeout = 1000000)
  870. })
  872. updated_neuronal_genes_ls_nowakowski_DivRG<- lapply(neuronal_genes_ls[12], function(x){
  873. UpdateSymbolList(x,
  874. timeout = 100000)
  875. })
  876. ## Astro
  877. updated_neuronal_genes_ls_nowakowski_astro<- lapply(neuronal_genes_ls[10], function(x){
  878. UpdateSymbolList(x,
  879. timeout = 20000)
  880. })
  882. ##OPC
  883. updated_neuronal_genes_ls_nowakowski_opc<- lapply(neuronal_genes_ls[11], function(x){
  884. UpdateSymbolList(x,
  885. timeout = 20000)
  886. })
  888. ## Other RG
  889. updated_neuronal_genes_ls_nowakowski_otherRG<- lapply(neuronal_genes_ls[13:16], function(x){
  890. UpdateSymbolList(x,
  891. timeout = 20000)
  892. })
  894. ## Enodthelial, pericytes etc.
  895. updated_neuronal_genes_ls_nowakowski_Endothelial<- lapply(neuronal_genes_ls[17:19], function(x){
  896. UpdateSymbolList(x,
  897. timeout = 20000)
  898. })
  900. ## MGE_newborn
  901. updated_neuronal_genes_ls_nowakowski_MGE_newborn<- lapply(neuronal_genes_ls[20:24], function(x){
  902. UpdateSymbolList(x,
  903. timeout = 20000)
  904. })
  906. ## IPCs
  907. updated_neuronal_genes_ls_nowakowski_IPC<- lapply(neuronal_genes_ls[25:29], function(x){
  908. UpdateSymbolList(x,
  909. timeout = 20000)
  910. })
  912. ## iNs
  913. updated_neuronal_genes_ls_nowakowski_iN<- lapply(neuronal_genes_ls[30:33], function(x){
  914. UpdateSymbolList(x,
  915. timeout = 20000)
  916. })
  918. ## eNs
  919. updated_neuronal_genes_ls_nowakowski_eN <- lapply(neuronal_genes_ls[c(34:43)], function(x){
  920. UpdateSymbolList(x,
  921. timeout = 20000)
  922. })
  924. load(file = paste0(datadirectory, "/reference_sheets/updated_neuronal_genes_ls_2.RData"))
  927. updated_neuronal_genes_ls <- c(updated_neuronal_genes_ls_liddelow, updated_neuronal_genes_ls_nowakowski_astro, updated_neuronal_genes_ls_nowakowski_Endothelial, updated_neuronal_genes_ls_nowakowski_iN, updated_neuronal_genes_ls_nowakowski_IPC, updated_neuronal_genes_ls_nowakowski_MGE_newborn, updated_neuronal_genes_ls_nowakowski_MGE_progenitors, updated_neuronal_genes_ls_nowakowski_eN, updated_neuronal_genes_ls_2)
  929. ## Couldn't update these gene sets
  930. updated_neuronal_genes_ls[["Div_RG_Sphase"]] <- neuronal_genes_ls[["Div_RG_Sphase"]]
  931. updated_neuronal_genes_ls[["Div_RG"]] <- neuronal_genes_ls[["Div_RG"]]
  933. names(updated_neuronal_genes_ls)[which(names(updated_neuronal_genes_ls) == "Astocyte")] <- "Astrocyte"
  935. updated_neuronal_genes_ls <- updated_neuronal_genes_ls[order(match(neuronal_genes_ls, updated_neuronal_genes_ls))]
  938. save(updated_neuronal_genes_ls, file = paste0(datadirectory, "/reference_sheets/updated_neuronal_genes_ls.RData"))
  939. ```

Neuronal_Glial_Types_Analysis_60D.Rmd at commit 1e1b492, no license · at the source

Overview

Authors: Zarina Greenberg1,2, Ella McDonald1,2, Alejandra Noreña Puerta1, Manam Inushi De Silva1,2, Cade Christensen1, Robert Adams1,2, Jenne Tran1,2, Paris Mazzachi1,2, Sebastian Loskarn1, Siti N Mubarokah2, Leanne Winner2, Drew Neavin3,4, Megan Maack5,6, Kristina L Elvidge5,6, Lisa Melton5, Mark R Hutchinson7,8, Kim M Hemsley2, Nicholas Smith9,10, Cedric Bardy1,2,11
ORCID iDs: Ella McDonald, Alejandra Noreña Puerta, Manam Inushi De Silva, Cade Christensen, Robert Adams, Paris Mazzachi, Sebastian Loskarn, Siti N Mubarokah, Leanne Winner, Drew Neavin, Megan Maack, Kristina L Elvidge, Lisa Melton, Mark R Hutchinson, Nicholas Smith, Cedric Bardy
  1. Laboratory for Human Neurophysiology and Genetics, South Australian Health and Medical Research Institute (SAHMRI), Adelaide, SA Australia
  2. Flinders Health and Medical Research Institute, College of Medicine and Public Health, Flinders University, Adelaide, SA Australia
  3. Translational Genomics, Garvan Institute of Medical Research, Darlinghurst, NSW Australia
  4. Centre for Population and Disease Genomics, Institute for Molecular Bioscience, The University of Queensland, Brisbane, QLD Australia
  5. Sanfilippo Children’s Foundation, Sydney, NSW Australia
  6. Childhood Dementia Initiative, Sydney, NSW Australia
  7. School of Pharmacy and Biomedical Science, Adelaide University, Adelaide, SA Australia
  8. Institute for Photonics and Advanced Sensing, Adelaide University, Adelaide, SA Australia
  9. Department of Neurology and Clinical Neurophysiology, Women’s and Children’s Health Network, Adelaide, SA Australia
  10. Paediatric Neurodegenerative Disease Research Group, Discipline of Paediatrics, College of Health, Adelaide University, Adelaide, SA Australia
  11. Brain Organoid Therapeutics, Adelaide, SA Australia
Institutions: Flinders University (Australia); South Australian Health and Medical Research Institute (Australia); Garvan Institute of Medical Research (Australia); The University of Queensland (Australia); Institute for Molecular Bioscience (Australia); Sanfilippo Children's Foundation (Australia); Adelaide University (Australia); The University of Adelaide (Australia); Women's and Children's Health Network (Australia)
Journal: Nature communications, volume 17, issue 1, article 9980
Dates: received 12 September 2024; accepted 5 August 2026; published online 20 August 2026
Type: Research article · Language: English
License: CC BY-NC-ND
Identifiers: DOI 10.1038/s41467-026-76837-1 · PMID 42764280 · PMCID PMC13590601 · OpenAlex W7203774183
Open access: gold, a free copy (OpenAlex)
Status: code verified
Categories: human (organism), other condition (population), Alzheimer's / dementia (population)
Methods: Spectral & time-frequency, Statistics, Smoothing, state filtering, decompositions, Machine learning, Preprocessing, Evoked potentials, Single-unit activity, calcium imaging, Physiology & signal measures
Keywords: Cellular neuroscience, Phenotypic screening, Induced pluripotent stem cells, Neurodegenerative diseases, Paediatric neurological disorders
MeSH: Dementia*, Machine Learning*, Mucopolysaccharidosis III*, Neuroprotective Agents*, Child, Drug Evaluation, Preclinical, Drug Repositioning, Female, Heparan Sulfate, Humans, Induced Pluripotent Stem Cells, Lysosomes (* major topic)
Topic: Lysosomal Storage Disorders Research (Physiology, Medicine), according to OpenAlex
Funding: Department of Health | National Health and Medical Research Council (NHMRC) (EPCD000025, MRF2024419)
Citations: not cited yet (Europe PMC); 100 references in the paper

Abstract

The abstract is not reproduced here: the paper's license (CC BY-NC-ND) does not allow it. Read it in the paper, at the publisher or on Europe PMC.

Repositories

Its files are read in the Code ↔ Paper reader above, with 23 matches between paragraphs and lines of code.

jacobgil/pytorch-grad-cam

License: MIT
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: 704393448a7b0c620ee2c2b9597723f1d7f17b3d, 13 August 2026
Languages: Python (54), Jupyter (7)
Size: 135 files, 61 scripts
Software Heritage: archived
Found in: the text, “Machine learning image models”
Holds: README, license file, environment (pyproject.toml, requirements.txt, setup.cfg, setup.py), tests, continuous integration, 7 notebooks
Not found: CITATION.cff, documentation
Tools: NumPy (41 files), PyTorch (37 files), OpenCV (18 files), Pillow (8 files), scikit-learn (2 files), SciPy (2 files), Matplotlib (1 file), Hugging Face Transformers (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
60 files
At the source: github.com/jacobgil/pytorch-grad-cam

bardylab/MPSIIIA_snRNAseq_DrugScreen_Paper_2026

License: none: the authors keep all their rights
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: 1e1b49286e5cf3122cb3a5cb44b2f9737be4ea02, 22 July 2026
Languages: R (6)
Size: 19 files, 6 scripts
Software Heritage: not archived
Found in: “Code availability”
Holds: README, 6 notebooks
Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: ggplot2 (6 files), Seurat (6 files), tidyverse (6 files), clusterProfiler (2 files), ComplexHeatmap (2 files), data.table (2 files), Harmony (2 files)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
7 files
At the source: github.com/bardylab/MPSIIIA_snRNAseq_DrugScreen_Paper_2026.git

bardylab/MPSIIIA_Machine_Learning_DrugScreen_Paper_2026

License: none: the authors keep all their rights
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: 2ba0e823f7d859fbef8e44829c8048a8b6e4c517, 20 August 2025
Languages: Python (14), Jupyter (2)
Size: 28 files, 16 scripts
Software Heritage: not archived
Found in: “Code availability”
Holds: README, 2 notebooks
Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: PyTorch (8 files), NumPy (7 files), pandas (6 files), Matplotlib (5 files), scikit-learn (5 files), seaborn (4 files), Pillow (3 files), XGBoost (3 files), SHAP (2 files), OpenCV (1 file), scikit-image (1 file), SciPy (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
17 files
At the source: github.com/bardylab/MPSIIIA_Machine_Learning_DrugScreen_Paper_2026

bardylab/MPSIIIA_snRNA-seq_DrugScreen_Paper_2026

License: none: the authors keep all their rights
State: the link is dead, verified on 27 September 2026
Evidence: found in the paper
Software Heritage: not archived
Found in: “Code availability”
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Availability: 1 check, the latest on 27 September 2026: the link is dead
  • 27 September 2026: the link is dead
At the source: github.com/bardylab/MPSIIIA_snRNA-seq_DrugScreen_Paper_2026.git

Code availability statement

The paper has a code availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

Read it in the paper: doi.org/10.1038/s41467-026-76837-1.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 4 repositories of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 80 scripts, each with its path and the digest of its content;
  • 23 matches between paragraphs of the paper and lines of the code (method lexical-v1);
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

No dataset and no data link were found in the paper.

Code and data availability statement

The paper has a code and data availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

Read it in the paper: doi.org/10.1038/s41467-026-76837-1.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 19 authors, 5 keywords, 12 MeSH terms, 1 funder, 99 references.

Cite

This paper

Greenberg, Z., McDonald, E., Noreña Puerta, A., De Silva, M. I., Christensen, C., Adams, R., Tran, J., Mazzachi, P., Loskarn, S., Mubarokah, S. N., Winner, L., Neavin, D., Maack, M., Elvidge, K. L., Melton, L., Hutchinson, M. R., Hemsley, K. M., Smith, N., & Bardy, C. (2026). Drug screen and machine learning predict neuroprotective agents in a preclinical human model of childhood dementia. Nature communications, 17(1), 9980. https://doi.org/10.1038/s41467-026-76837-1

BibTeX

@article{greenberg2026drug,
author = {Greenberg, Zarina and McDonald, Ella and Noreña Puerta, Alejandra and De Silva, Manam Inushi and Christensen, Cade and Adams, Robert and Tran, Jenne and Mazzachi, Paris and Loskarn, Sebastian and Mubarokah, Siti N and Winner, Leanne and Neavin, Drew and Maack, Megan and Elvidge, Kristina L and Melton, Lisa and Hutchinson, Mark R and Hemsley, Kim M and Smith, Nicholas and Bardy, Cedric},
title = {{Drug screen and machine learning predict neuroprotective agents in a preclinical human model of childhood dementia}},
journal = {Nature communications},
year = {2026},
month = aug,
volume = {17},
number = {1},
pages = {9980},
publisher = {Nature Publishing Group},
issn = {2041-1723},
doi = {10.1038/s41467-026-76837-1},
url = {https://doi.org/10.1038/s41467-026-76837-1},
pmid = {42764280},
pmcid = {PMC13590601}
}

RIS

TY - JOUR
AU - Greenberg, Zarina
AU - McDonald, Ella
AU - Noreña Puerta, Alejandra
AU - De Silva, Manam Inushi
AU - Christensen, Cade
AU - Adams, Robert
AU - Tran, Jenne
AU - Mazzachi, Paris
AU - Loskarn, Sebastian
AU - Mubarokah, Siti N
AU - Winner, Leanne
AU - Neavin, Drew
AU - Maack, Megan
AU - Elvidge, Kristina L
AU - Melton, Lisa
AU - Hutchinson, Mark R
AU - Hemsley, Kim M
AU - Smith, Nicholas
AU - Bardy, Cedric
TI - Drug screen and machine learning predict neuroprotective agents in a preclinical human model of childhood dementia
T2 - Nature communications
J2 - Nat Commun
PY - 2026
DA - 2026/08/20
VL - 17
IS - 1
SP - 9980
SN - 2041-1723
PB - Nature Publishing Group
DO - 10.1038/s41467-026-76837-1
UR - https://doi.org/10.1038/s41467-026-76837-1
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s41467-026-76837-1",
"type": "article-journal",
"title": "Drug screen and machine learning predict neuroprotective agents in a preclinical human model of childhood dementia",
"container-title": "Nature communications",
"author": [
{
"family": "Greenberg",
"given": "Zarina"
},
{
"family": "McDonald",
"given": "Ella"
},
{
"family": "Noreña Puerta",
"given": "Alejandra"
},
{
"family": "De Silva",
"given": "Manam Inushi"
},
{
"family": "Christensen",
"given": "Cade"
},
{
"family": "Adams",
"given": "Robert"
},
{
"family": "Tran",
"given": "Jenne"
},
{
"family": "Mazzachi",
"given": "Paris"
},
{
"family": "Loskarn",
"given": "Sebastian"
},
{
"family": "Mubarokah",
"given": "Siti N"
},
{
"family": "Winner",
"given": "Leanne"
},
{
"family": "Neavin",
"given": "Drew"
},
{
"family": "Maack",
"given": "Megan"
},
{
"family": "Elvidge",
"given": "Kristina L"
},
{
"family": "Melton",
"given": "Lisa"
},
{
"family": "Hutchinson",
"given": "Mark R"
},
{
"family": "Hemsley",
"given": "Kim M"
},
{
"family": "Smith",
"given": "Nicholas"
},
{
"family": "Bardy",
"given": "Cedric"
}
],
"container-title-short": "Nat Commun",
"volume": "17",
"issue": "1",
"page": "9980",
"DOI": "10.1038/s41467-026-76837-1",
"PMID": "42764280",
"PMCID": "PMC13590601",
"ISSN": "2041-1723",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s41467-026-76837-1",
"language": "en",
"issued": {
"date-parts": [
[
2026,
8,
20
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1038/s41467-026-71112-9 [code]
Modelling synaptic dysfunction in childhood dementia using human iPSC-derived cortical networks.
Journal: Nature communications
In common: Harmony, clusterProfiler, Seurat, 3 other tools, Alzheimer's / dementia, 34 references, 13 authors
[2] doi:10.1038/s42003-026-10957-8 [code]
Brain defence by the extracellular matrix protein Cochlin.
Journal: Communications biology
In common: SHAP, Hugging Face Transformers, XGBoost, 16 other tools
[3] doi:10.1016/j.xcrm.2026.102766 [code]
A longitudinal single-cell and spatial multiomic atlas of pediatric high-grade glioma.
Journal: Cell reports. Medicine
In common: Harmony, clusterProfiler, ComplexHeatmap, 13 other tools, other condition, 1 reference
[4] doi:10.1002/imt2.70163 [code]
Spatial multi-omics unveils sphingolipid metabolic reprogramming within the retinal pathological niche.
Journal: iMeta
In common: Harmony, clusterProfiler, ComplexHeatmap, 13 other tools
[5] doi:10.1186/s13073-026-01704-z [code]
Gene expression profiling enables refined parcellation of cortical layers in the heterogeneous human cerebral cortex.
Journal: Genome medicine
In common: Harmony, XGBoost, ComplexHeatmap, 12 other tools
[6] doi:10.1038/s41586-026-10214-2 [code]
Multidimensional profiling of heterogeneity in supratentorial ependymomas.
Journal: Nature
In common: Harmony, clusterProfiler, ComplexHeatmap, 10 other tools, other condition, 1 reference
[7] doi:10.1038/s41467-026-71595-6 [code]
A single-cell and spatial atlas of early human olfactory development.
Journal: Nature communications
In common: Harmony, ComplexHeatmap, Seurat, 6 other tools, 4 references
[8] doi:10.1038/s41514-026-00391-9 [code]
Region-specific transcriptional signatures of brain aging in the absence of neuropathology at the single-cell level.
Journal: npj aging
In common: Harmony, ComplexHeatmap, Seurat, 10 other tools, 1 reference
[9] doi:10.1016/j.xcrm.2026.102651 [code]
Integrative CSF profiling identifies disease-specific immune responses in leptomeningeal disease.
Journal: Cell reports. Medicine
In common: Harmony, ComplexHeatmap, Seurat, 10 other tools, other condition
[10] doi:10.1126/sciadv.aeg3223 [code]
The extreme diversity of retinal amacrine cells has deep evolutionary roots.
Journal: Science advances
In common: XGBoost, clusterProfiler, ComplexHeatmap, 10 other tools

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.