rnaSentry operates on a SummarizedExperiment: the
expression matrix in the assay slot and the sample metadata in
colData. Most wet-lab and clinical users, however, start
with two plain text files:
This vignette shows the exact bridge between those two files and the
pipeline: read the files, keep only the samples present in both (the
single most common pitfall), validate the intake with
load_counts(), audit the design with
design_audit(), and run the full pipeline.
The example below writes two small synthetic files to a temporary
directory as stand-ins for your counts.csv and
clinical.csv; replace the file paths and column names with
your own and the block works unchanged.
library(rnaSentry)
library(SummarizedExperiment)
set.seed(42)
n_genes <- 300
n_samples <- 30
genes <- paste0("GENE", seq_len(n_genes))
samples <- paste0("S", seq_len(n_samples))
# Raw integer counts (bulk RNA-seq style), genes x samples.
cnt <- matrix(rpois(n_genes * n_samples, lambda = 200),
nrow = n_genes, ncol = n_samples,
dimnames = list(genes, samples))
# Clinical metadata, one row per sample.
clin <- S4Vectors::DataFrame(
time = round(rexp(n_samples, rate = 0.08), 1),
event = rbinom(n_samples, size = 1, prob = 0.5),
batch = factor(rep(c("B1", "B2", "B3"), length.out = n_samples)),
stage = factor(sample(c("I", "II", "III"), n_samples, replace = TRUE)),
row.names = samples
)
# --- Deliberately introduce the two classic file problems -------------------
# A sample present in the counts file but not the clinical file,
cnt <- cbind(cnt, S31 = cnt[, "S2"])
# and a clinical row with no matching counts column.
clin <- rbind(clin, S32 = clin["S1", ])
# Plus two intake issues load_counts() should catch: a non-integer count...
cnt["GENE1", "S1"] <- cnt["GENE1", "S1"] + 0.5
# ...and a gene with no counts anywhere.
cnt["GENE300", ] <- 0L
td <- tempfile("csvdemo")
dir.create(td)
write.csv(cnt, file.path(td, "counts.csv"))
write.csv(clin, file.path(td, "clinical.csv"))
list.files(td)
#> [1] "clinical.csv" "counts.csv"Every reader sees the same two files a real user exports from their
pipeline: counts.csv (genes in rows, samples in columns)
and clinical.csv (one row per sample).
counts <- as.matrix(read.csv(file.path(td, "counts.csv"), row.names = 1))
clinical <- read.csv(file.path(td, "clinical.csv"), row.names = 1)
dim(counts)
#> [1] 300 31
dim(clinical)
#> [1] 31 4Note the as.matrix(): load_counts() expects
a matrix, and read.csv() returns a
data.frame.
The most common pitfall is sample IDs that do not match between the
two files: one file may contain samples the other lacks (a failed RNA
extraction, a sample re-labeled, or a version mismatch). Passing them to
load_counts() un-checked produces a hard error or, worse,
silently mis-aligned columns. The fix is a single
intersect() that keeps only samples found in
both files:
n_counts_before <- ncol(counts)
n_clin_before <- nrow(clinical)
common <- intersect(colnames(counts), rownames(clinical))
counts <- counts[, common, drop = FALSE]
clinical <- clinical[common, , drop = FALSE]
data.frame(
"counts columns" = c(before = n_counts_before, after = ncol(counts)),
"clinical rows" = c(before = n_clin_before, after = nrow(clinical)),
row.names = c("before intersect", "after intersect")
)
#> counts.columns clinical.rows
#> before intersect 31 31
#> after intersect 30 30This example planted one sample in each file that the other lacks, so
the 30 truly shared samples survive and
S31/S32 are dropped. Always check the
dimensions after intersecting: if you lose hundreds of samples, suspect
a version or delimiter mismatch between the two files.
load_counts() constructs the object and records every
intake issue it finds in a flag ledger stored in
metadata(se)$flags:
se <- load_counts(counts, clinical)
#> Warning: Count matrix contains non-integer values. Downstream stages expect
#> integer counts.
#> Warning: 1 gene(s) have all-zero or all-NA counts and will be retained but may
#> be filtered later.
se
#> class: SummarizedExperiment
#> dim: 300 30
#> metadata(1): flags
#> assays(1): counts
#> rownames(300): GENE1 GENE2 ... GENE299 GENE300
#> rowData names(0):
#> colnames(30): S1 S2 ... S29 S30
#> colData names(4): time event batch stageS4Vectors::metadata(se)$flags
#> check severity detail
#> 1 non_integer_counts warning Count matrix contains non-integer values.
#> 2 empty_rows warning 1 gene(s) have all-zero or all-NA counts.
#> stage
#> 1 load_counts
#> 2 load_countsThe ledger flags the two problems we planted (the non-integer count and the all-zero gene) and carries them into the rendered report’s audit trail. Duplicate sample IDs, missing metadata, and Ensembl-style gene IDs are caught the same way. The pipeline still runs on a flagged object – the ledger’s job is to make sure you saw the issue first.
design_audit() tests every candidate covariate against
the leading expression principal component and flags potential
confounders, before any per-gene modeling is invested in:
da <- design_audit(se, design_vars = c("batch", "stage"))
#> Warning in stats::chisq.test(tab): Chi-squared approximation may be incorrect
#> Warning in stats::chisq.test(tab): Chi-squared approximation may be incorrect
da
#> rnaSentry design audit: 2 candidate variable(s), surrogate PC1.
#> Recommended design formula: ~ batch + stage
#> 0 variable(s) flagged as potential confounders.
#> Interaction LRT run; 0 variable(s) with significant surrogate interaction.
#> No issues flagged.run_rnaSentry() orchestrates design audit, PCA audit,
signature construction with repeated cross-validation, and the
survival-modeling stages on the same object. It needs a numeric
follow-up time column and a 0/1
event indicator in colData:
run <- run_rnaSentry(se, time_col = "time", event_col = "event",
top_n = 5, repeats = 1, folds = 2, seed = 42,
design_vars = c("batch", "stage"),
render_report = FALSE)
#> Warning in stats::chisq.test(tab): Chi-squared approximation may be incorrect
#> Warning in stats::chisq.test(tab): Chi-squared approximation may be incorrect
#> Warning: Only 14 event(s) for a 5-gene signature (2.8 events per parameter),
#> below 'min_events_per_parameter' = 5. Cross-validated concordance and hazard
#> ratios are unstable at this event count; consider fewer genes or a larger
#> cohort.
run
#> rnaSentry pipeline run.
#> - design_audit
#> - pca_audit
#> - build_signature
#> - km_curve
#> - cox_model
#> - survival_parametricOn 30 samples the events_per_parameter guardrail warns
that cross-validated concordance is unstable – exactly as it should:
this cohort is a demonstration, not a discovery set. See
vignette("case-study-small-cohort") for what that warning
means and ?generate_report for how to render the audited
HTML report of a run.
The run locked the discovered signature: run_rnaSentry()
is single-use per session, so a second call in the same session is
refused until the lock is released. This is the reproducibility
guardrail against silently re-selecting signature genes after survival
analysis.
counts vs logcounts
(exprs for microarray)load_counts() stores your matrix under the assay name
counts, and the pipeline computes
log2(counts + 1) when it runs. If your input file is
already normalized (CPM, TPM, VST), store it under the name
logcounts instead so it is used as-is and not
double-logged. For microarray log2 intensities, use the conventional
assay name exprs, which the pipeline likewise uses
as-is:
Either way the bridge is identical; only the assay name tells the pipeline which scale your data is in.
sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04.1 LTS
#>
#> Matrix products: default
#> BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so; LAPACK version 3.12.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: Etc/UTC
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats4 stats graphics grDevices utils datasets methods
#> [8] base
#>
#> other attached packages:
#> [1] survival_3.8-12 SummarizedExperiment_1.43.0
#> [3] Biobase_2.73.2 GenomicRanges_1.65.4
#> [5] Seqinfo_1.3.2 IRanges_2.47.5
#> [7] S4Vectors_0.51.10 BiocGenerics_0.59.12
#> [9] generics_0.1.4 MatrixGenerics_1.25.0
#> [11] matrixStats_1.5.0 rnaSentry_0.99.19
#> [13] BiocStyle_2.41.0
#>
#> loaded via a namespace (and not attached):
#> [1] gtable_0.3.6 ggplot2_4.0.3 xfun_0.61
#> [4] bslib_0.12.0 httr2_1.3.0 lattice_0.23-1
#> [7] tzdb_0.5.0 vctrs_0.7.3 tools_4.6.1
#> [10] curl_8.0.0 tibble_3.3.1 R.oo_1.27.1
#> [13] pkgconfig_2.0.3 Matrix_1.7-6 data.table_1.18.6.1
#> [16] RColorBrewer_1.1-3 rentrez_1.2.4 S7_0.2.2
#> [19] lifecycle_1.0.5 farver_2.1.2 compiler_4.6.1
#> [22] statmod_1.5.2 htmltools_0.5.9 sys_3.4.3
#> [25] buildtools_1.0.0 sass_0.4.10 yaml_2.3.12
#> [28] pillar_1.11.1 jquerylib_0.1.4 tidyr_1.3.2
#> [31] R.utils_2.13.0 DelayedArray_0.39.8 cachem_1.1.0
#> [34] limma_3.99.0 abind_1.4-8 tidyselect_1.2.1
#> [37] digest_0.6.39 dplyr_1.2.1 purrr_1.2.2
#> [40] labeling_0.4.3 maketools_1.3.2 splines_4.6.1
#> [43] fastmap_1.2.0 grid_4.6.1 GEOquery_2.81.28
#> [46] cli_3.6.6 SparseArray_1.13.4 magrittr_2.0.5
#> [49] S4Arrays_1.13.2 XML_3.99-0.25 readr_2.2.0
#> [52] withr_3.0.3 scales_1.4.0 rmarkdown_2.32
#> [55] XVector_0.53.0 otel_0.2.0 R.methodsS3_1.8.2
#> [58] hms_1.1.4 evaluate_1.0.5 knitr_1.52
#> [61] rlang_1.3.0 glue_1.8.1 BiocManager_1.30.27
#> [64] xml2_1.6.0 jsonlite_2.0.0 R6_2.6.1