---
title: "2. Getting started: from two CSVs to a SummarizedExperiment"
author: "Abdelkarim Hani Ghrieb"
date: "`r Sys.Date()`"
output:
  BiocStyle::html_document:
    toc: true
    toc_float: true
vignette: >
  %\VignetteIndexEntry{2. Getting started: from two CSVs to a SummarizedExperiment}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r getting-started-setup, include = FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
```

# Motivation

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:

- a **counts matrix** (genes in rows, samples in columns), and
- a **clinical table** (one row per sample, with follow-up time and event
  plus any batch, stage, or treatment covariates).

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.

# Build two example CSV files

```{r build-csv}
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)
```

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).

# Step 1: read both files

```{r read-csv}
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)
dim(clinical)
```

Note the `as.matrix()`: `load_counts()` expects a matrix, and
`read.csv()` returns a `data.frame`.

# Step 2: intersect the sample IDs

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:

```{r intersect-samples}
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")
)
```

This 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.

# Step 3: build and validate the SummarizedExperiment

`load_counts()` constructs the object and records every intake issue it
finds in a flag ledger stored in `metadata(se)$flags`:

```{r load-counts}
se <- load_counts(counts, clinical)
se
```

```{r load-flags}
S4Vectors::metadata(se)$flags
```

The 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.

# Step 4: audit the design before modeling

`design_audit()` tests every candidate covariate against the leading
expression principal component and flags potential confounders, before any
per-gene modeling is invested in:

```{r design-audit}
da <- design_audit(se, design_vars = c("batch", "stage"))
da
```

# Step 5: run the full pipeline

`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`:

```{r run-pipeline, warning = TRUE}
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)
run
```

On 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.

```{r unlock-signature}
run$stages$build_signature <-
    lock_signature(run$stages$build_signature, lock = FALSE)
```

# Assay naming: `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:

```{r assay-naming}
se_norm <- se
SummarizedExperiment::assayNames(se_norm) <- "logcounts"
```

Either way the bridge is identical; only the assay name tells the pipeline
which scale your data is in.

# Session info

```{r session-info}
sessionInfo()
```