Skip to content

Commit

Permalink
Merge branch 'main' of https://github.com/jr-leary7/scLANE into dev
Browse files Browse the repository at this point in the history
  • Loading branch information
jr-leary7 committed Nov 13, 2023
2 parents 6fc5a5a + 4766369 commit 19acd5f
Show file tree
Hide file tree
Showing 7 changed files with 116 additions and 21 deletions.
137 changes: 116 additions & 21 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,7 @@
<!-- badges: start -->

[![R-CMD-check](https://github.com/jr-leary7/scLANE/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/jr-leary7/scLANE/actions/workflows/R-CMD-check.yaml)
[![Bioc-check](https://github.com/jr-leary7/scLANE/actions/workflows/bioc-check.yaml/badge.svg)](https://github.com/jr-leary7/scLANE/actions/workflows/bioc-check.yaml)
![release](https://img.shields.io/github/v/release/jr-leary7/scLANE?color=purple)
![last
commit](https://img.shields.io/github/last-commit/jr-leary7/scLANE/main?color=darkgreen)
Expand Down Expand Up @@ -94,11 +95,8 @@ conflict.

``` r
library(dplyr)
library(scater)
library(scLANE)
library(ggplot2)
select <- dplyr::select
filter <- dplyr::filter
```

## Input data
Expand All @@ -117,8 +115,101 @@ The PCA embeddings show us a pretty simple trajectory that’s strongly
correlated with the first principal component.

``` r
plotPCA(sim_data, colour_by = "cell_time_normed") +
data.frame(sim_data@int_colData$reducedDims@listData$PCA[, 1:2]) %>%
mutate(pseudotime = sim_data$cell_time_normed) %>%
ggplot(aes(x = PC1, y = PC2, color = pseudotime)) +
geom_point(size = 2, alpha = 0.75, stroke = 0) +
scale_color_gradientn(colors = viridisLite::plasma(n = 20)) +
labs(x = "PC 1", y = "PC 2", color = "Pseudotime") +
theme_scLANE(umap = TRUE)
#> Loading required package: S4Vectors
#> Loading required package: stats4
#> Loading required package: BiocGenerics
#>
#> Attaching package: 'BiocGenerics'
#> The following objects are masked from 'package:dplyr':
#>
#> combine, intersect, setdiff, union
#> The following objects are masked from 'package:stats':
#>
#> IQR, mad, sd, var, xtabs
#> The following objects are masked from 'package:base':
#>
#> anyDuplicated, aperm, append, as.data.frame, basename, cbind,
#> colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
#> get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
#> match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
#> Position, rank, rbind, Reduce, rownames, sapply, setdiff, sort,
#> table, tapply, union, unique, unsplit, which.max, which.min
#>
#> Attaching package: 'S4Vectors'
#> The following objects are masked from 'package:dplyr':
#>
#> first, rename
#> The following object is masked from 'package:utils':
#>
#> findMatches
#> The following objects are masked from 'package:base':
#>
#> expand.grid, I, unname
#> Loading required package: SingleCellExperiment
#> Loading required package: SummarizedExperiment
#> Loading required package: MatrixGenerics
#> Loading required package: matrixStats
#>
#> Attaching package: 'matrixStats'
#> The following object is masked from 'package:dplyr':
#>
#> count
#>
#> Attaching package: 'MatrixGenerics'
#> The following objects are masked from 'package:matrixStats':
#>
#> colAlls, colAnyNAs, colAnys, colAvgsPerRowSet, colCollapse,
#> colCounts, colCummaxs, colCummins, colCumprods, colCumsums,
#> colDiffs, colIQRDiffs, colIQRs, colLogSumExps, colMadDiffs,
#> colMads, colMaxs, colMeans2, colMedians, colMins, colOrderStats,
#> colProds, colQuantiles, colRanges, colRanks, colSdDiffs, colSds,
#> colSums2, colTabulates, colVarDiffs, colVars, colWeightedMads,
#> colWeightedMeans, colWeightedMedians, colWeightedSds,
#> colWeightedVars, rowAlls, rowAnyNAs, rowAnys, rowAvgsPerColSet,
#> rowCollapse, rowCounts, rowCummaxs, rowCummins, rowCumprods,
#> rowCumsums, rowDiffs, rowIQRDiffs, rowIQRs, rowLogSumExps,
#> rowMadDiffs, rowMads, rowMaxs, rowMeans2, rowMedians, rowMins,
#> rowOrderStats, rowProds, rowQuantiles, rowRanges, rowRanks,
#> rowSdDiffs, rowSds, rowSums2, rowTabulates, rowVarDiffs, rowVars,
#> rowWeightedMads, rowWeightedMeans, rowWeightedMedians,
#> rowWeightedSds, rowWeightedVars
#> Loading required package: GenomicRanges
#> Loading required package: IRanges
#>
#> Attaching package: 'IRanges'
#> The following objects are masked from 'package:dplyr':
#>
#> collapse, desc, slice
#> Loading required package: GenomeInfoDb
#>
#> Attaching package: 'GenomicRanges'
#> The following object is masked from 'package:magrittr':
#>
#> subtract
#> Loading required package: Biobase
#> Welcome to Bioconductor
#>
#> Vignettes contain introductory material; view with
#> 'browseVignettes()'. To cite Bioconductor, see
#> 'citation("Biobase")', and for packages 'citation("pkgname")'.
#>
#> Attaching package: 'Biobase'
#> The following object is masked from 'package:MatrixGenerics':
#>
#> rowMedians
#> The following objects are masked from 'package:matrixStats':
#>
#> anyMissing, rowMedians
#> Registered S3 method overwritten by 'SparseArray':
#> method from
#> rowsum.dgCMatrix DelayedArray
```

<img src="man/figures/README-plot-sims-pt-1.png" width="100%" />
Expand All @@ -127,7 +218,11 @@ We also see that the data are not clustered by subject, which indicates
that gene dynamics are mostly homogeneous across subjects.

``` r
plotPCA(sim_data, colour_by = "subject") +
data.frame(sim_data@int_colData$reducedDims@listData$PCA[, 1:2]) %>%
mutate(subject = sim_data$subject) %>%
ggplot(aes(x = PC1, y = PC2, color = subject)) +
geom_point(size = 2, alpha = 0.75, stroke = 0) +
labs(x = "PC 1", y = "PC 2", color = "Subject ID") +
theme_scLANE(umap = TRUE)
```

Expand Down Expand Up @@ -181,7 +276,7 @@ scLANE_models_glm <- testDynamic(sim_data,
#> Registered S3 method overwritten by 'bit':
#> method from
#> print.ri gamlss
#> scLANE testing completed for 100 genes across 1 lineage in 35.114 secs
#> scLANE testing completed for 100 genes across 1 lineage in 29.977 secs
```

After the function finishes running, we use `getResultsDE()` to generate
Expand All @@ -200,13 +295,13 @@ select(scLANE_res_glm, Gene, Lineage, Test_Stat, P_Val, P_Val_Adj, Gene_Dynamic_
col.names = c("Gene", "Lineage", "LRT stat.", "P-value", "Adj. p-value", "Predicted dynamic status"))
```

| Gene | Lineage | LRT stat. | P-value | Adj. p-value | Predicted dynamic status |
|:---------|:--------|----------:|--------:|-------------:|-------------------------:|
| MFSD2B | A | 209.755 | 0.000 | 0.000 | 1 |
| GOLGA8EP | A | 4.390 | 0.036 | 0.537 | 0 |
| NOP14 | A | 7.712 | 0.005 | 0.132 | 0 |
| TMCO3 | A | 167.582 | 0.000 | 0.000 | 1 |
| BATF2 | A | 5.216 | 0.074 | 0.537 | 0 |
| Gene | Lineage | LRT stat. | P-value | Adj. p-value | Predicted dynamic status |
|:-------|:--------|----------:|--------:|-------------:|-------------------------:|
| MFSD2B | A | 209.755 | 0.000 | 0.000 | 1 |
| SMG1 | A | 5.062 | 0.080 | 0.419 | 0 |
| UAP1L1 | A | 9.118 | 0.010 | 0.199 | 0 |
| TMCO3 | A | 167.582 | 0.000 | 0.000 | 1 |
| MYOF | A | 4.015 | 0.045 | 0.419 | 0 |

### GEE mode

Expand All @@ -229,7 +324,7 @@ scLANE_models_gee <- testDynamic(sim_data,
id.vec = sim_data$subject,
cor.structure = "ar1",
n.cores = 4)
#> scLANE testing completed for 100 genes across 1 lineage in 2.201 mins
#> scLANE testing completed for 100 genes across 1 lineage in 3.076 mins
```

We again generate the table of DE test results. The variance of the
Expand All @@ -248,11 +343,11 @@ select(scLANE_res_gee, Gene, Lineage, Test_Stat, P_Val, P_Val_Adj, Gene_Dynamic_

| Gene | Lineage | Wald stat. | P-value | Adj. p-value | Predicted dynamic status |
|:---------|:--------|-----------:|--------:|-------------:|-------------------------:|
| DGUOK | A | 200675.460 | 0.000 | 0.000 | 1 |
| VDAC1 | A | 13.025 | 0.001 | 0.025 | 0 |
| CKAP4 | A | 48922.644 | 0 | 0 | 1 |
| TBCC | A | 32.151 | 0 | 0 | 1 |
| GOLGA8EP | A | NA | NA | NA | 0 |
| WAPAL | A | 3254.351 | 0.000 | 0.000 | 1 |
| IDH3G | A | 863.479 | 0.000 | 0.000 | 1 |
| PFDN2 | A | 2131.763 | 0 | 0 | 1 |
| IDH3G | A | 863.479 | 0 | 0 | 1 |

### GLMM mode

Expand All @@ -276,7 +371,7 @@ scLANE_models_glmm <- testDynamic(sim_data,
glmm.adaptive = TRUE,
id.vec = sim_data$subject,
n.cores = 4)
#> scLANE testing completed for 100 genes across 1 lineage in 3.427 mins
#> scLANE testing completed for 100 genes across 1 lineage in 4.863 mins
```

{% note %}
Expand Down Expand Up @@ -306,8 +401,8 @@ select(scLANE_res_glmm, Gene, Lineage, Test_Stat, P_Val, P_Val_Adj, Gene_Dynamic

| Gene | Lineage | LRT stat. | P-value | Adj. p-value | Predicted dynamic status |
|:--------|:--------|----------:|--------:|-------------:|-------------------------:|
| DDX1 | A | 131.707 | 0.000 | 0 | 1 |
| TSPAN1 | A | 78.616 | 0.000 | 0 | 1 |
| PGLS | A | 129.086 | 0.000 | 0 | 1 |
| TSPAN1 | A | 85.987 | 0.000 | 0 | 1 |
| WDSUB1 | A | NA | NA | NA | 0 |
| FAM135B | A | NA | NA | NA | 0 |
| NLGN4Y | A | 9.878 | 0.627 | 1 | 0 |
Expand Down
Binary file modified man/figures/README-plot-models-gee-1.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified man/figures/README-plot-models-glm-1.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified man/figures/README-plot-models-glmm-1.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified man/figures/README-plot-sims-pt-1.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified man/figures/README-plot-sims-subj-1.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified man/figures/README-unnamed-chunk-2-1.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.

0 comments on commit 19acd5f

Please sign in to comment.