> **Status:** `draft`<br>
> **Template class:** `SOURCE-BACKED WORKFLOW`

## Purpose

Prepare a named list of modality-appropriate numeric matrices for MOFA2. The
source workflow used a strict intersection of samples across views and selected
high-variance features. This page keeps that transparent baseline while making
the number and names of views user-controlled.

## Inputs and parameters

The input `.qs2` file must contain a named list. Each element is a matrix with
features in rows and biologically meaningful sample IDs in columns. The matrices
must already be prepared at a modality-appropriate scale, or the explicit
per-view transformation must be added beside the corresponding input before
the checks below. This page does not guess transformations from view names.

```{r}
#| label: parameters

input_views_path <- "input/modality_views.qs2"
prepared_views_path <- "output/mofa/prepared_views.qs2"
view_manifest_path <- "output/mofa/view_manifest.tsv"
sample_overlap_path <- "output/mofa/sample_overlap.tsv"
feature_manifest_path <- "output/mofa/feature_manifest.tsv"

# Empty means keep every eligible non-constant feature in every view.
# Add only the views for which a top-N rule is scientifically justified.
# Example: top_features <- list(transcriptome = 5000L, proteome = NA_integer_)
top_features <- list()
```

::: {.callout-note}
## View-specific preprocessing

The source workflow used logCPM for count-derived RNA/miRNA views, while
methylation, protein and binary views were prepared according to their own
measurement scales. Those choices are examples, not rules based on a view name.
For example, a count-derived view may be prepared explicitly with
`edgeR::DGEList()`, `edgeR::calcNormFactors(method = "TMMwsp")` and
`edgeR::cpm(log = TRUE, prior.count = 1)`. Do not apply bulk-RNA logCPM to a
non-count modality, automatically z-score every view, or turn missing values
into zero. Consider the data distribution, intended likelihood, feature
variance and technical artifacts before constructing the list.
:::

## Validate orientation and identifiers

```{r}
#| label: validate-view-contract

dir.create(dirname(prepared_views_path), recursive = TRUE, showWarnings = FALSE)
views <- qs2::qs_read(input_views_path)

if (!is.list(views) || is.null(names(views)) || any(!nzchar(names(views)))) {
  stop("The input must be a non-empty named list of feature-by-sample matrices.")
}
if (anyDuplicated(names(views))) {
  stop("View names must be unique.")
}
if (length(views) < 2L) {
  stop("A multi-omics MOFA workflow requires at least two views.")
}

sample_sets <- lapply(views, function(x) {
  if (!is.matrix(x) || !is.numeric(x)) {
    stop("Every view must be a numeric matrix; do not silently coerce it.")
  }
  if (nrow(x) == 0L || ncol(x) == 0L) {
    stop("Every view must have at least one feature and one sample.")
  }
  if (is.null(rownames(x)) || any(!nzchar(rownames(x))) ||
      anyDuplicated(rownames(x))) {
    stop("Every view needs non-empty, unique feature names.")
  }
  if (is.null(colnames(x)) || any(!nzchar(colnames(x))) ||
      anyDuplicated(colnames(x))) {
    stop("Every view needs non-empty, unique sample IDs.")
  }
  if (any(!is.na(x) & !is.finite(x))) {
    stop("Views may contain NA values, but not Inf or NaN values.")
  }
  if (any(rowSums(!is.na(x)) == 0L)) {
    stop("Remove or explicitly handle all-missing features before MOFA.")
  }
  colnames(x)
})

view_names <- names(views)
all_sample_ids <- sort(unique(unlist(sample_sets, use.names = FALSE)))
sample_overlap <- do.call(rbind, lapply(view_names, function(view_name) {
  data.frame(
    view = view_name,
    sample_id = all_sample_ids,
    present = all_sample_ids %in% sample_sets[[view_name]],
    stringsAsFactors = FALSE
  )
}))

common_samples <- Reduce(intersect, sample_sets)
if (length(common_samples) < 2L) {
  stop("Fewer than two samples are present in every selected view.")
}

view_manifest <- do.call(rbind, lapply(view_names, function(view_name) {
  present <- sample_sets[[view_name]]
  other_samples <- unlist(sample_sets[setdiff(view_names, view_name)],
                          use.names = FALSE)
  data.frame(
    view = view_name,
    samples_before = length(present),
    samples_retained = sum(present %in% common_samples),
    samples_dropped = sum(!present %in% common_samples),
    unique_to_view = paste(setdiff(present, other_samples), collapse = ","),
    stringsAsFactors = FALSE
  )
}))

readr::write_tsv(view_manifest, view_manifest_path)
readr::write_tsv(sample_overlap, sample_overlap_path)
```

::: {.callout-note}
## Sample matching

The source workflow uses samples present in every selected view as the simple
canonical baseline. This is transparent and easy to reason about, but it is not
the only possible MOFA design. MOFA2 can model missing values rather than using
a hidden imputation step. For deliberately incomplete multi-omics designs,
consult the current MOFA2 documentation and encode missingness explicitly
instead of automatically imputing values or silently dropping samples.
:::

## Match samples and select features

```{r}
#| label: match-and-filter-views

feature_manifest_rows <- vector("list", length(view_names))
for (view_name in view_names) {
  x <- views[[view_name]][, common_samples, drop = FALSE]
  feature_variance <- apply(x, 1L, stats::var, na.rm = TRUE)
  eligible_features <- rownames(x)[is.finite(feature_variance) &
                                   feature_variance > 0]

  requested <- if (view_name %in% names(top_features)) {
    top_features[[view_name]]
  } else {
    NA_integer_
  }
  if (length(requested) != 1L ||
      (!is.na(requested) &&
       (requested < 1L || requested != as.integer(requested)))) {
    stop("Each top_features value must be one positive integer or NA.")
  }

  if (is.na(requested)) {
    selected_features <- eligible_features
  } else {
    ordered_features <- eligible_features[
      order(feature_variance[eligible_features], decreasing = TRUE)
    ]
    selected_features <- head(ordered_features, requested)
  }
  if (length(selected_features) == 0L) {
    stop("A view has no eligible non-constant features after filtering.")
  }

  views[[view_name]] <- x[selected_features, , drop = FALSE]
  feature_manifest_rows[[view_name]] <- data.frame(
    view = view_name,
    features_before = nrow(x),
    features_after = nrow(views[[view_name]]),
    requested_top_features = requested,
    retained_nonconstant_features = length(eligible_features),
    internal_na_values = sum(is.na(x)),
    stringsAsFactors = FALSE
  )
}
feature_manifest <- do.call(rbind, feature_manifest_rows)

if (!all(vapply(views, function(x) identical(colnames(x), common_samples),
                logical(1)))) {
  stop("All prepared views must use the same matched sample order.")
}
if (any(vapply(views, function(x) anyDuplicated(rownames(x)) > 0L,
                 logical(1)))) {
  stop("Duplicate feature names remain after filtering.")
}

readr::write_tsv(feature_manifest, feature_manifest_path)
qs2::qs_save(object = views, file = prepared_views_path)

print(view_manifest)
print(feature_manifest)
cat("Prepared", length(views), "views and", length(common_samples),
    "matched samples.\n")
```

## Output contract

The saved object is a named list of numeric feature-by-sample matrices with a
common, explicitly checked sample order. The manifests account for sample
matching and feature filtering. Internal `NA` values are preserved; no
imputation or NA-to-zero conversion is performed here.
