1  Skin: QC & Filtering

1.1 Set up Seurat workspace

# Load libraries
library(data.table)
library(devtools)
library(presto)
library(glmGamPoi)
library(sctransform)
library(Seurat)
library(tidyverse)
library(miQC)   
library(SeuratWrappers)
library(flexmix)
library(SingleCellExperiment)
library(SummarizedExperiment)
library(readxl)
library(fishpond)
library(Matrix)
library(speckle)
library(scater)
library(patchwork)
library(vctrs)
library(alevinQC)
library(harmony)
library(scDblFinder)
library(cellXY)

# Set global options for Seurat v5 objects
options(Seurat.object.assay.version = 'v5')

1.2 Load alevin-fry output

Note using U+S+A counts via snRNA

P101_Skin_Post3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P101_Skin_Post3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 80763 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P101_Skin_Pre3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P101_Skin_Pre3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 76442 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P103_Skin_Post3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P103_Skin_Post3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 64399 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P103_Skin_Pre3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P103_Skin_Pre3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 59420 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P104_Skin_Post3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P104_Skin_Post3rd_2.5mgIpi_RNA_fixed020224_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 57061 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P104_Skin_Pre3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P104_Skin_Pre3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 63314 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P105_Skin_Post3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P105_Skin_Post3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 65472 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P105_Skin_Pre3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P105_Skin_Pre3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 70513 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P106_Skin_Post3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P106_Skin_Post3rd_2.5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 65665 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P106_Skin_Pre3rd_2.5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P106_Skin_Pre3rd_2.5mgIpi_RNA_fixed021324_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 65342 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P108_Skin_Post3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P108_Skin_Post3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 69381 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P108_Skin_Pre3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P108_Skin_Pre3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 69324 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P109_Skin_Post3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P109_Skin_Post3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 71971 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P109_Skin_Pre3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P109_Skin_Pre3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 75117 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P110_Skin_Post3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P110_Skin_Post3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 66843 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P110_Skin_Pre3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P110_Skin_Pre3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 45085 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P111_Skin_Post3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P111_Skin_Post3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 63149 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

P111_Skin_Pre3rd_5mgIpi_RNA <- loadFry(fryDir = "/jsimonlab/projects/Wu/Melanoma_scRNA_Eryn/alevinfry_072023/P111_Skin_Pre3rd_5mgIpi_RNA_alevin_quant_crlikeem/", outputFormat = "snRNA")

locating quant file

Reading meta data

USA mode: TRUE

Processing 62703 genes and 67303 barcodes

Using pre-defined output format: snrna

Building the 'counts' assay, which contains U S A

Constructing output SingleCellExperiment object

Done

1.3 Convert gene IDs to gene symbols, collapse transcripts to gene symbols

tx2gene <- read.table("/jsimonlab/genomes/hg38/gencode.v43.annotation_splici/gene_id_to_name.tsv",header=F,sep="\t",col.names=c("tx","gene"))

# Use rownames of first object to apply to all
exp.txId <- rownames(P101_Skin_Post3rd_2.5mgIpi_RNA)
exp.geneId <- as.vector(tx2gene$gene[match(exp.txId, tx2gene$tx)])
exp.tx.grp <- t(sparse.model.matrix(~ 0 + exp.geneId))


P101_Skin_Post3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P101_Skin_Post3rd_2.5mgIpi_RNA)
rownames(P101_Skin_Post3rd_2.5mgIpi_RNA.summarized) <- rownames(P101_Skin_Post3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P101_Skin_Pre3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P101_Skin_Pre3rd_2.5mgIpi_RNA)
rownames(P101_Skin_Pre3rd_2.5mgIpi_RNA.summarized) <- rownames(P101_Skin_Pre3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P103_Skin_Post3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P103_Skin_Post3rd_2.5mgIpi_RNA)
rownames(P103_Skin_Post3rd_2.5mgIpi_RNA.summarized) <- rownames(P103_Skin_Post3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P103_Skin_Pre3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P103_Skin_Pre3rd_2.5mgIpi_RNA)
rownames(P103_Skin_Pre3rd_2.5mgIpi_RNA.summarized) <- rownames(P103_Skin_Pre3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P104_Skin_Post3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P104_Skin_Post3rd_2.5mgIpi_RNA)
rownames(P104_Skin_Post3rd_2.5mgIpi_RNA.summarized) <- rownames(P104_Skin_Post3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P104_Skin_Pre3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P104_Skin_Pre3rd_2.5mgIpi_RNA)
rownames(P104_Skin_Pre3rd_2.5mgIpi_RNA.summarized) <- rownames(P104_Skin_Pre3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P105_Skin_Post3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P105_Skin_Post3rd_2.5mgIpi_RNA)
rownames(P105_Skin_Post3rd_2.5mgIpi_RNA.summarized) <- rownames(P105_Skin_Post3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P105_Skin_Pre3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P105_Skin_Pre3rd_2.5mgIpi_RNA)
rownames(P105_Skin_Pre3rd_2.5mgIpi_RNA.summarized) <- rownames(P105_Skin_Pre3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P106_Skin_Post3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P106_Skin_Post3rd_2.5mgIpi_RNA)
rownames(P106_Skin_Post3rd_2.5mgIpi_RNA.summarized) <- rownames(P106_Skin_Post3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P106_Skin_Pre3rd_2.5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P106_Skin_Pre3rd_2.5mgIpi_RNA)
rownames(P106_Skin_Pre3rd_2.5mgIpi_RNA.summarized) <- rownames(P106_Skin_Pre3rd_2.5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P108_Skin_Post3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P108_Skin_Post3rd_5mgIpi_RNA)
rownames(P108_Skin_Post3rd_5mgIpi_RNA.summarized) <- rownames(P108_Skin_Post3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P108_Skin_Pre3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P108_Skin_Pre3rd_5mgIpi_RNA)
rownames(P108_Skin_Pre3rd_5mgIpi_RNA.summarized) <- rownames(P108_Skin_Pre3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P109_Skin_Post3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P109_Skin_Post3rd_5mgIpi_RNA)
rownames(P109_Skin_Post3rd_5mgIpi_RNA.summarized) <- rownames(P109_Skin_Post3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P109_Skin_Pre3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P109_Skin_Pre3rd_5mgIpi_RNA)
rownames(P109_Skin_Pre3rd_5mgIpi_RNA.summarized) <- rownames(P109_Skin_Pre3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P110_Skin_Post3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P110_Skin_Post3rd_5mgIpi_RNA)
rownames(P110_Skin_Post3rd_5mgIpi_RNA.summarized) <- rownames(P110_Skin_Post3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P110_Skin_Pre3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P110_Skin_Pre3rd_5mgIpi_RNA)
rownames(P110_Skin_Pre3rd_5mgIpi_RNA.summarized) <- rownames(P110_Skin_Pre3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P111_Skin_Post3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P111_Skin_Post3rd_5mgIpi_RNA)
rownames(P111_Skin_Post3rd_5mgIpi_RNA.summarized) <- rownames(P111_Skin_Post3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

P111_Skin_Pre3rd_5mgIpi_RNA.summarized <- exp.tx.grp %*% counts(P111_Skin_Pre3rd_5mgIpi_RNA)
rownames(P111_Skin_Pre3rd_5mgIpi_RNA.summarized) <- rownames(P111_Skin_Pre3rd_5mgIpi_RNA.summarized) %>% str_replace_all(".+.geneId","")

1.4 Create Seurat objects

P101_Skin_Post3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P101_Skin_Post3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P101_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P101_Skin_Pre3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P103_Skin_Post3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P103_Skin_Post3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P103_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P103_Skin_Pre3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P104_Skin_Post3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P104_Skin_Post3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P104_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P104_Skin_Pre3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P105_Skin_Post3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P105_Skin_Post3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P105_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P105_Skin_Pre3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P106_Skin_Post3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P106_Skin_Post3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P106_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- CreateSeuratObject(P106_Skin_Pre3rd_2.5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P108_Skin_Post3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P108_Skin_Post3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P108_Skin_Pre3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P108_Skin_Pre3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P109_Skin_Post3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P109_Skin_Post3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P109_Skin_Pre3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P109_Skin_Pre3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P110_Skin_Post3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P110_Skin_Post3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P110_Skin_Pre3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P110_Skin_Pre3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P111_Skin_Post3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P111_Skin_Post3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

P111_Skin_Pre3rd_5mgIpi_RNA.seurat <- CreateSeuratObject(P111_Skin_Pre3rd_5mgIpi_RNA.summarized)

Warning: Feature names cannot have underscores ('_'), replacing with dashes
('-')

1.5 Add sample name prefix to CBs

P101_Skin_Post3rd_2.5mgIpi_RNA.seurat <- RenameCells(P101_Skin_Post3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P101_Skin_Post3rd_2.5mgIpi_RNA")
P101_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- RenameCells(P101_Skin_Pre3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P101_Skin_Pre3rd_2.5mgIpi_RNA")
P103_Skin_Post3rd_2.5mgIpi_RNA.seurat <- RenameCells(P103_Skin_Post3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P103_Skin_Post3rd_2.5mgIpi_RNA")
P103_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- RenameCells(P103_Skin_Pre3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P103_Skin_Pre3rd_2.5mgIpi_RNA")
P104_Skin_Post3rd_2.5mgIpi_RNA.seurat <- RenameCells(P104_Skin_Post3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P104_Skin_Post3rd_2.5mgIpi_RNA")
P104_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- RenameCells(P104_Skin_Pre3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P104_Skin_Pre3rd_2.5mgIpi_RNA")
P105_Skin_Post3rd_2.5mgIpi_RNA.seurat <- RenameCells(P105_Skin_Post3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P105_Skin_Post3rd_2.5mgIpi_RNA")
P105_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- RenameCells(P105_Skin_Pre3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P105_Skin_Pre3rd_2.5mgIpi_RNA")
P106_Skin_Post3rd_2.5mgIpi_RNA.seurat <- RenameCells(P106_Skin_Post3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P106_Skin_Post3rd_2.5mgIpi_RNA")
P106_Skin_Pre3rd_2.5mgIpi_RNA.seurat <- RenameCells(P106_Skin_Pre3rd_2.5mgIpi_RNA.seurat,add.cell.id = "P106_Skin_Pre3rd_2.5mgIpi_RNA")
P108_Skin_Post3rd_5mgIpi_RNA.seurat <- RenameCells(P108_Skin_Post3rd_5mgIpi_RNA.seurat,add.cell.id = "P108_Skin_Post3rd_5mgIpi_RNA")
P108_Skin_Pre3rd_5mgIpi_RNA.seurat <- RenameCells(P108_Skin_Pre3rd_5mgIpi_RNA.seurat,add.cell.id = "P108_Skin_Pre3rd_5mgIpi_RNA")
P109_Skin_Post3rd_5mgIpi_RNA.seurat <- RenameCells(P109_Skin_Post3rd_5mgIpi_RNA.seurat,add.cell.id = "P109_Skin_Post3rd_5mgIpi_RNA")
P109_Skin_Pre3rd_5mgIpi_RNA.seurat <- RenameCells(P109_Skin_Pre3rd_5mgIpi_RNA.seurat,add.cell.id = "P109_Skin_Pre3rd_5mgIpi_RNA")
P110_Skin_Post3rd_5mgIpi_RNA.seurat <- RenameCells(P110_Skin_Post3rd_5mgIpi_RNA.seurat,add.cell.id = "P110_Skin_Post3rd_5mgIpi_RNA")
P110_Skin_Pre3rd_5mgIpi_RNA.seurat <- RenameCells(P110_Skin_Pre3rd_5mgIpi_RNA.seurat,add.cell.id = "P110_Skin_Pre3rd_5mgIpi_RNA")
P111_Skin_Post3rd_5mgIpi_RNA.seurat <- RenameCells(P111_Skin_Post3rd_5mgIpi_RNA.seurat,add.cell.id = "P111_Skin_Post3rd_5mgIpi_RNA")
P111_Skin_Pre3rd_5mgIpi_RNA.seurat <- RenameCells(P111_Skin_Pre3rd_5mgIpi_RNA.seurat,add.cell.id = "P111_Skin_Pre3rd_5mgIpi_RNA")

1.6 Merge objects

P109_Skin_Post3rd was removed due to low viable cell count

merged.18279.skin <- merge(x = P101_Skin_Pre3rd_2.5mgIpi_RNA.seurat, y=c(P101_Skin_Post3rd_2.5mgIpi_RNA.seurat, P103_Skin_Pre3rd_2.5mgIpi_RNA.seurat, P103_Skin_Post3rd_2.5mgIpi_RNA.seurat, P104_Skin_Post3rd_2.5mgIpi_RNA.seurat, P104_Skin_Pre3rd_2.5mgIpi_RNA.seurat, P105_Skin_Post3rd_2.5mgIpi_RNA.seurat, P105_Skin_Pre3rd_2.5mgIpi_RNA.seurat, P106_Skin_Pre3rd_2.5mgIpi_RNA.seurat, P106_Skin_Post3rd_2.5mgIpi_RNA.seurat, P108_Skin_Pre3rd_5mgIpi_RNA.seurat, P108_Skin_Post3rd_5mgIpi_RNA.seurat, P109_Skin_Pre3rd_5mgIpi_RNA.seurat, P110_Skin_Pre3rd_5mgIpi_RNA.seurat, P110_Skin_Post3rd_5mgIpi_RNA.seurat, P111_Skin_Pre3rd_5mgIpi_RNA.seurat, P111_Skin_Post3rd_5mgIpi_RNA.seurat))
dim(merged.18279.skin)

[1]   61217 1124593

1.7 QC filter

Use relatively loose initial filters

merged.18279.skin <- subset(merged.18279.skin, subset = nCount_RNA > 500 & nFeature_RNA > 250)
dim(merged.18279.skin)

[1] 61217 64133

merged.18279.skin <- PercentageFeatureSet(merged.18279.skin, pattern = "^MT-", col.name = "percent.mt")
merged.18279.skin <- RunMiQC(merged.18279.skin, 
                        percent.mt = "percent.mt", 
                        nFeature_RNA = "nFeature_RNA", 
                        posterior.cutoff = 0.6, 
                        model.slot = "flexmix_model")

Warning: Adding a command log without an assay associated with it

merged.18279.skin <- subset(merged.18279.skin, miQC.keep == "keep")

dim(merged.18279.skin)

[1] 61217 46891

data.frame(table(str_replace_all(colnames(merged.18279.skin),"_RNA_.+","")))

                         Var1 Freq
1  P101_Skin_Post3rd_2.5mgIpi 2338
2   P101_Skin_Pre3rd_2.5mgIpi  622
3  P103_Skin_Post3rd_2.5mgIpi 5340
4   P103_Skin_Pre3rd_2.5mgIpi 1214
5  P104_Skin_Post3rd_2.5mgIpi 7396
6   P104_Skin_Pre3rd_2.5mgIpi 4291
7  P105_Skin_Post3rd_2.5mgIpi 4668
8   P105_Skin_Pre3rd_2.5mgIpi 3341
9  P106_Skin_Post3rd_2.5mgIpi 2214
10  P106_Skin_Pre3rd_2.5mgIpi 2223
11   P108_Skin_Post3rd_5mgIpi 3601
12    P108_Skin_Pre3rd_5mgIpi 1277
13    P109_Skin_Pre3rd_5mgIpi 2144
14   P110_Skin_Post3rd_5mgIpi 1979
15    P110_Skin_Pre3rd_5mgIpi 1593
16   P111_Skin_Post3rd_5mgIpi 1949
17    P111_Skin_Pre3rd_5mgIpi  701

1.8 Add meta data

patient <- str_split_i(colnames(merged.18279.skin),"_",1)
site <- str_split_i(colnames(merged.18279.skin),"_",2)
timepoint <- str_split_i(colnames(merged.18279.skin),"_",3)
IpiCohort <- str_split_i(colnames(merged.18279.skin),"_",4)
assay <- str_split_i(colnames(merged.18279.skin),"_",5)
barcode <- str_split_i(colnames(merged.18279.skin),"_",6)
sample <- paste0(patient,"_",site,"_",timepoint,"_",IpiCohort,"_",assay)

merged.18279.skin <- AddMetaData(merged.18279.skin, patient, col.name="Patient")
merged.18279.skin <- AddMetaData(merged.18279.skin, site, col.name="Site")
merged.18279.skin <- AddMetaData(merged.18279.skin, timepoint, col.name="Timepoint")
merged.18279.skin <- AddMetaData(merged.18279.skin, IpiCohort, col.name="IpiCohort")
merged.18279.skin <- AddMetaData(merged.18279.skin, assay, col.name="Assay")
merged.18279.skin <- AddMetaData(merged.18279.skin, sample, col.name="Sample")

1.9 Join then re-split RNA counts layers by Sample

merged.18279.skin[['RNA']] <- JoinLayers(merged.18279.skin[['RNA']])
merged.18279.skin[["RNA"]] <- split(merged.18279.skin[["RNA"]], f = merged.18279.skin$Sample)

1.10 Save merged object

saveRDS(merged.18279.skin,"Skin_scRNA_Part1.rds")

1.11 Get session info

sessionInfo()

R version 4.3.1 (2023-06-16)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Rocky Linux 8.10 (Green Obsidian)

Matrix products: default
BLAS/LAPACK: /usr/lib64/libopenblasp-r0.3.15.so;  LAPACK version 3.9.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: America/New_York
tzcode source: system (glibc)

attached base packages:
[1] stats4    stats     graphics  grDevices utils     datasets  methods  
[8] base     

other attached packages:
 [1] cellXY_0.99.0               scDblFinder_1.14.0         
 [3] harmony_1.2.0               alevinQC_1.16.1            
 [5] vctrs_0.6.5                 patchwork_1.3.0            
 [7] scater_1.28.0               scuttle_1.10.3             
 [9] speckle_1.0.0               Matrix_1.6-4               
[11] fishpond_2.6.2              readxl_1.4.3               
[13] SingleCellExperiment_1.22.0 SummarizedExperiment_1.30.2
[15] Biobase_2.60.0              GenomicRanges_1.52.1       
[17] GenomeInfoDb_1.36.4         IRanges_2.34.1             
[19] S4Vectors_0.38.2            BiocGenerics_0.46.0        
[21] MatrixGenerics_1.12.3       matrixStats_1.2.0          
[23] flexmix_2.3-19              lattice_0.22-5             
[25] SeuratWrappers_0.3.19       miQC_1.8.0                 
[27] lubridate_1.9.3             forcats_1.0.0              
[29] stringr_1.5.1               dplyr_1.1.4                
[31] purrr_1.0.2                 readr_2.1.5                
[33] tidyr_1.3.1                 tibble_3.2.1               
[35] ggplot2_3.4.4               tidyverse_2.0.0            
[37] Seurat_5.1.0                SeuratObject_5.0.2         
[39] sp_2.1-3                    sctransform_0.4.1          
[41] glmGamPoi_1.12.2            presto_1.0.0               
[43] Rcpp_1.0.12                 devtools_2.4.5             
[45] usethis_2.2.2               data.table_1.15.0          

loaded via a namespace (and not attached):
  [1] fs_1.6.3                  spatstat.sparse_3.0-3    
  [3] bitops_1.0-7              httr_1.4.7               
  [5] RColorBrewer_1.1-3        profvis_0.3.8            
  [7] tools_4.3.1               utf8_1.2.4               
  [9] R6_2.5.1                  DT_0.31                  
 [11] lazyeval_0.2.2            uwot_0.1.16              
 [13] urlchecker_1.0.1          withr_3.0.0              
 [15] GGally_2.2.1              gridExtra_2.3            
 [17] progressr_0.14.0          cli_3.6.2                
 [19] spatstat.explore_3.2-6    fastDummies_1.7.3        
 [21] spatstat.data_3.0-4       ggridges_0.5.6           
 [23] pbapply_1.7-2             Rsamtools_2.16.0         
 [25] R.utils_2.12.3            parallelly_1.37.0        
 [27] sessioninfo_1.2.2         limma_3.56.2             
 [29] RSQLite_2.3.5             BiocIO_1.10.0            
 [31] generics_0.1.3            gtools_3.9.5             
 [33] ica_1.0-3                 spatstat.random_3.2-2    
 [35] ggbeeswarm_0.7.2          fansi_1.0.6              
 [37] abind_1.4-5               R.methodsS3_1.8.2        
 [39] lifecycle_1.0.4           yaml_2.3.8               
 [41] edgeR_3.42.4              Rtsne_0.17               
 [43] blob_1.2.4                grid_4.3.1               
 [45] dqrng_0.3.2               promises_1.2.1           
 [47] crayon_1.5.2              shinydashboard_0.7.2     
 [49] miniUI_0.1.1.1            beachmat_2.16.0          
 [51] cowplot_1.1.3             KEGGREST_1.40.1          
 [53] metapod_1.8.0             pillar_1.9.0             
 [55] knitr_1.45                rjson_0.2.21             
 [57] xgboost_1.7.7.1           future.apply_1.11.1      
 [59] codetools_0.2-19          leiden_0.4.3.1           
 [61] glue_1.7.0                remotes_2.4.2.1          
 [63] png_0.1-8                 spam_2.10-0              
 [65] org.Mm.eg.db_3.18.0       cellranger_1.1.0         
 [67] gtable_0.3.4              cachem_1.0.8             
 [69] xfun_0.42                 S4Arrays_1.2.0           
 [71] mime_0.12                 survival_3.5-8           
 [73] statmod_1.5.0             bluster_1.10.0           
 [75] ellipsis_0.3.2            fitdistrplus_1.1-11      
 [77] ROCR_1.0-11               nlme_3.1-164             
 [79] bit64_4.0.5               RcppAnnoy_0.0.22         
 [81] irlba_2.3.5.1             vipor_0.4.7              
 [83] KernSmooth_2.23-22        DBI_1.2.2                
 [85] colorspace_2.1-0          nnet_7.3-19              
 [87] tidyselect_1.2.0          bit_4.0.5                
 [89] compiler_4.3.1            BiocNeighbors_1.18.0     
 [91] DelayedArray_0.26.7       plotly_4.10.4            
 [93] rtracklayer_1.60.1        scales_1.3.0             
 [95] lmtest_0.9-40             digest_0.6.34            
 [97] goftest_1.2-3             spatstat.utils_3.0-4     
 [99] rmarkdown_2.25            XVector_0.40.0           
[101] htmltools_0.5.7           pkgconfig_2.0.3          
[103] sparseMatrixStats_1.12.2  fastmap_1.1.1            
[105] rlang_1.1.3               htmlwidgets_1.6.4        
[107] shiny_1.8.0               DelayedMatrixStats_1.22.6
[109] farver_2.1.1              zoo_1.8-12               
[111] jsonlite_1.8.8            BiocParallel_1.34.2      
[113] R.oo_1.26.0               BiocSingular_1.16.0      
[115] RCurl_1.98-1.14           magrittr_2.0.3           
[117] modeltools_0.2-23         GenomeInfoDbData_1.2.10  
[119] dotCall64_1.1-1           munsell_0.5.0            
[121] viridis_0.6.5             reticulate_1.35.0        
[123] stringi_1.8.3             zlibbioc_1.46.0          
[125] MASS_7.3-60.0.1           org.Hs.eg.db_3.18.0      
[127] plyr_1.8.9                pkgbuild_1.4.3           
[129] ggstats_0.5.1             parallel_4.3.1           
[131] listenv_0.9.1             ggrepel_0.9.5            
[133] deldir_2.0-2              Biostrings_2.68.1        
[135] splines_4.3.1             tensor_1.5               
[137] hms_1.1.3                 locfit_1.5-9.8           
[139] igraph_2.0.2              spatstat.geom_3.2-8      
[141] RcppHNSW_0.6.0            reshape2_1.4.4           
[143] ScaledMatrix_1.8.1        pkgload_1.3.4            
[145] XML_3.99-0.16.1           evaluate_0.23            
[147] scran_1.28.2              BiocManager_1.30.22      
[149] tzdb_0.4.0                httpuv_1.6.14            
[151] RANN_2.6.1                polyclip_1.10-6          
[153] future_1.33.1             scattermore_1.2          
[155] rsvd_1.0.5                xtable_1.8-4             
[157] restfulr_0.0.15           svMisc_1.2.3             
[159] RSpectra_0.16-1           later_1.3.2              
[161] viridisLite_0.4.2         AnnotationDbi_1.64.1     
[163] GenomicAlignments_1.36.0  memoise_2.0.1            
[165] beeswarm_0.4.0            tximport_1.28.0          
[167] cluster_2.1.6             timechange_0.3.0         
[169] globals_0.16.2