Part 1 — Linear Methods for Multi-Omic Integration

This workshop notebook compares two simple linear strategies for integrating TCGA-BRCA multi-omics data:

1. Early integration: concatenate transcriptomics, proteomics, and methylation features, then fit one linear classifier.

2. Late integration: fit one linear classifier per omic, then aggregate their predictions.

The data are already processed, aligned, and complete. The patient ID is the index, and the supervised target is the subtype column.

Learning objectives

By the end of this notebook you should be able to:

• Load already processed omics matrices.

• Use a single shared train/test split across all integration methods.

• Visualize how omics differ using PCA.

• Train single-omic linear baselines.

• Train and evaluate early integration by concatenation.

• Train and evaluate late integration by prediction averaging.

• Discuss why simple integration motivates patient-level multi-omic representations in the next part.

1. Import Libraries and Helpers

# ── Imports ──────────────────────────────────────────────────────────────────
from pathlib import Path
import pickle

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

from IPython.display import display

from sklearn.base import clone
from sklearn.decomposition import PCA
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import (
    accuracy_score,
    balanced_accuracy_score,
    classification_report,
    confusion_matrix,
)
from sklearn.model_selection import train_test_split
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import LabelEncoder, StandardScaler

# ── Custom Imports ────────────────────────────────────────────────────────────
from s2_helpers import load_omics, evaluate_predictions, plot_confusion_matrix

# ── Reproducibility ───────────────────────────────────────────────────────────
RANDOM_STATE = 42

2. Load the prepared omics data and patient train / test splits

We separate:

• subtype as the prediction target

• all remaining columns as omic features

• the index as the patient ID

DATA_DIR = Path("/data/")

X_omics , y = load_omics(DATA_DIR , omic_keys=['transcriptomics' , 'proteomics', 'methylation'])

with open(f'{DATA_DIR}/patient_splits.pkl' , 'rb') as file : 
    data = pickle.load(file)

train_ids = data['train_ids']
test_ids = data['test_ids']

Omic view dimensions:
  transcriptomics:  500 patients ×  29995 features
  proteomics     :  500 patients ×    464 features
  methylation    :  500 patients × 200000 features

Subtype counts:

paper_BRCA_Subtype_PAM50
LumA      237
LumB      100
Basal      97
Her2       41
Normal     25
Name: count, dtype: int64

3. PCA view of each omic

PCA is not used for the classifiers below. It is only a quick visual check: if the panels separate patients differently, they may contain complementary information.

label_encoder = LabelEncoder()
y_encoded = label_encoder.fit_transform(y)
class_names = list(label_encoder.classes_)


def plot_pca_by_omic(X_omics, y_encoded, class_names):
    """Plot the first two PCs for each omic using all patients."""
    for name, X in X_omics.items():
        pca_pipe = Pipeline([
            ("scaler", StandardScaler()),
            ("pca", PCA(n_components=2, random_state=RANDOM_STATE)),
        ])
        scores = pca_pipe.fit_transform(X)
        evr = pca_pipe.named_steps["pca"].explained_variance_ratio_

        plt.figure(figsize=(5.5, 4.5))
        scatter = plt.scatter(scores[:, 0], scores[:, 1], c=y_encoded, s=25, alpha=0.8)
        plt.title(f"{name}: PCA")
        plt.xlabel(f"PC1 ({evr[0]:.1%} variance)")
        plt.ylabel(f"PC2 ({evr[1]:.1%} variance)")
        handles, _ = scatter.legend_elements()
        plt.legend(handles, class_names, title="subtype", bbox_to_anchor=(1.02, 1), loc="upper left")
        plt.tight_layout()
        plt.show()

plot_pca_by_omic(X_omics, y_encoded, class_names)

Interpretation (click to expand)

Transcriptomics shows the clearest visible subtype structure: Basal (purple) forms a distinct cluster separated from the LumA/LumB/Normal samples. This suggests that transcriptomic variation contains strong subtype-related signal, which may translate into stronger performance for transcriptomics-based models, although this must be confirmed using supervised prediction.

Proteomics shows weaker visible subtype structure in the first two principal components: samples from different subtypes largely overlap, with only a few outliers separated along PC1. This suggests that subtype-related variation is less dominant in the major variance directions, which may make classification more challenging, although predictive performance must be evaluated directly.

Methylation shows weaker partial structure: some separation along PC1 is visible (Basal/Her2 trending right, LumA/LumB trending left), but subtype groups remain substantially overlapping compared with transcriptomics.

Low variance explained everywhere (7–16% per PC) means these plots are only 2D projections of high-dimensional data. Fine cluster shapes should not be over-interpreted; only broad patterns, such as the stronger subtype separation in transcriptomics, could be considered meaningful.

4. Model helpers

We use the same regularized linear classifier everywhere. Standardization is included inside the pipeline so scaling is learned only from the training data.

def make_linear_classifier(): 
    """Regularized linear model for high-dimensional omics classification."""
    return Pipeline([
        ("scaler", StandardScaler()),
        ("model", LogisticRegression(
            C=1.0,
            l1_ratio=0,        # equivalent to L2
            solver="lbfgs", #instead of liblinear because it doesn't handle 3+ classes
            class_weight="balanced",
            max_iter=2000,
            random_state=RANDOM_STATE,
        ),
        ),
    ])

def fit_predict(model, X_train, y_train, X_test):
    """Fit a model and return class predictions and probabilities for the test set."""
    fitted = clone(model).fit(X_train, y_train)
    pred = pd.Series(fitted.predict(X_test), index=X_test.index, name="prediction")
    proba = pd.DataFrame(fitted.predict_proba(X_test), index=X_test.index, columns=fitted.classes_)
    return fitted, pred, proba

5. Single-omic benchmarks

Why integrate in the first place?

Before integrating data types, measure how far each omic gets on its own..

Before comparing integration strategies, we need a reason to integrate in the first place. Each omic layer below is trained as a standalone linear classifier, using the same train/test split, model, and metrics as every method that follows.

The question isn’t “which single omic is best?”; it’s whether the omics agree with each other or not. If every omic already predicted the same subtype for every patient, integration would add complexity without adding information. If they disagree on a meaningful fraction of patients, that disagreement is exactly the signal integration is trying to recover.

results = []
single_omic_models = {}
single_omic_probabilities = {}

base_model = make_linear_classifier()

for name, X in X_omics.items():
    model, pred, proba = fit_predict(
        base_model,
        X.loc[train_ids],
        y[train_ids],
        X.loc[test_ids],
    )
    single_omic_models[name] = model
    single_omic_probabilities[name] = proba
    results.append(evaluate_predictions(y[test_ids], pred, f"single omic: {name}"))

results_df = pd.DataFrame(results).sort_values("balanced_accuracy", ascending=False)
results_df

single omic: transcriptomics
────────────────────────────
  Accuracy          : 0.824 
  Balanced accuracy : 0.799

single omic: proteomics
───────────────────────
  Accuracy          : 0.720 
  Balanced accuracy : 0.643

single omic: methylation
────────────────────────
  Accuracy          : 0.712 
  Balanced accuracy : 0.688

	model	accuracy	balanced_accuracy
0	single omic: transcriptomics	0.824	0.798667
2	single omic: methylation	0.712	0.688333
1	single omic: proteomics	0.720	0.643333

# How often do the single-omic models actually agree on a prediction?
single_omic_preds = {
    name: single_omic_models[name].predict(X_omics[name].loc[test_ids])
    for name in X_omics
}
pred_table = pd.DataFrame(single_omic_preds, index=test_ids)
pred_table["true_label"] = y[test_ids]

n_unique_preds = pred_table[list(X_omics)].nunique(axis=1)
agreement_summary = pd.DataFrame({
    "n_patients": [
        (n_unique_preds == 1).sum(),
        (n_unique_preds > 1).sum(),
    ],
}, index=["all omics agree", "omics disagree"])
agreement_summary["pct"] = (100 * agreement_summary["n_patients"] / len(test_ids)).round(1)

print(f"Test patients where every omic predicts the same subtype: "
      f"{agreement_summary.loc['all omics agree', 'pct']}%")
display(agreement_summary)

# Look at a few disputed patients: does any single omic get the true label right
# where others miss it? This is the concrete case for integration.
disputed = pred_table.loc[n_unique_preds > 1].copy()
disputed["any_omic_correct"] = disputed[list(X_omics)].eq(disputed["true_label"], axis=0).any(axis=1)
print(f"Of the disputed patients, {disputed['any_omic_correct'].mean():.0%} "
      f"have at least one omic that got the true subtype right.")
disputed.head(10)

Test patients where every omic predicts the same subtype: 60.0%

	n_patients	pct
all omics agree	75	60.0
omics disagree	50	40.0

Of the disputed patients, 92% have at least one omic that got the true subtype right.

	transcriptomics	proteomics	methylation	true_label	any_omic_correct
TCGA-E9-A5UO	LumA	LumB	LumA	LumB	True
TCGA-D8-A1XY	LumB	LumA	LumB	LumA	True
TCGA-D8-A27G	LumA	LumA	LumB	LumA	True
TCGA-AR-A5QQ	Basal	Basal	Normal	Basal	True
TCGA-BH-A1FN	Basal	Basal	LumB	LumB	True
TCGA-E2-A108	Normal	LumA	LumA	Normal	True
TCGA-AO-A0JC	LumB	LumA	LumA	LumA	True
TCGA-AR-A0TX	Her2	LumA	Her2	Her2	True
TCGA-BH-A0BJ	Normal	LumA	LumA	LumA	True
TCGA-AR-A24R	LumB	LumA	LumB	LumB	True

Interpretation (click to expand)

We asked each omic to predict the subtype alone. On 65% of patients, all three omics agree. These represent cases where each modality provides consistent information, so integration may provide less additional benefit compared with cases where the omics disagree.

In the remaining 35%, the predictions disagree, suggesting that different omics capture different aspects of the underlying biology or contain different levels of uncertainty for those patients.

In 89% of disagreement cases, at least one omic-specific model correctly predicted the subtype.The signal was not absent; it was distributed across modalities, with no single omic capturing all of the relevant information consistently.

Therefore, the goal is not simply to identify the “best” omic, but to build models that can combine these partial views into an integrated representation that captures information from multiple sources.

7. Early integration by concatenation

What’s the simplest way to to combine multi omics, and at what does this simplicity cost us?

Early integration joins all omic features into one wide matrix and trains one classifier.

This is easy and often strong, but it increases dimensionality and does not explicitly distinguish shared, redundant, and omic-specific signal.

def prefix_columns(X: pd.DataFrame, omic_name: str) -> pd.DataFrame:
    """Add an omic prefix so features remain traceable after concatenation."""
    X_prefixed = X.copy()
    X_prefixed.columns = [f"{omic_name}::{col}" for col in X.columns]
    return X_prefixed

X_early = pd.concat(
    [prefix_columns(X, name) for name, X in X_omics.items()],
    axis=1,
)

print(f"Early integration matrix: {X_early.shape[0]} patients x {X_early.shape[1]} features")

early_model, early_pred, early_proba = fit_predict(
    base_model,
    X_early.loc[train_ids],
    y[train_ids],
    X_early.loc[test_ids],
)

results.append(evaluate_predictions(y[test_ids], early_pred, "early integration: concatenation"))
results_df = pd.DataFrame(results).sort_values("balanced_accuracy", ascending=False)
results_df

Early integration matrix: 500 patients x 230459 features

early integration: concatenation
────────────────────────────────
  Accuracy          : 0.752 
  Balanced accuracy : 0.702

	model	accuracy	balanced_accuracy
0	single omic: transcriptomics	0.824	0.798667
3	early integration: concatenation	0.752	0.702333
2	single omic: methylation	0.712	0.688333
1	single omic: proteomics	0.720	0.643333

Interpretation (click to expand)

Transcriptomics alone is the best performer (0.669 balanced accuracy), outperforming the other single-omic models and the concatenated model.

Early integration underperforms transcriptomics alone (0.647 vs 0.669), showing that simply combining all modalities does not guarantee improved prediction. Additional omic features may introduce redundant or noisy information that the model must separate from the strongest signal.

8. Inspect feature usage in the early-integration model

When early integration works, does it actually use all omics or lean on one?

Concatenating features gives the model access to everything at once, but access isn’t the same as use. Before trusting an early-integration result, it’s worth checking which omics its largest coefficients actually come from. A linear model on concatenated features can be inspected directly through coefficient magnitude, so this is a cheap, concrete check rather than a guess: does the top-100 coefficient list spread across omics, or does one omic quietly dominate?

This matters for what comes next: if early integration turns out to lean on a single omic, that is itself an argument for late integration (Section 9 onward) –> keeping omics as separate models makes this kind of imbalance visible and controllable (via weights) instead of hidden inside one shared coefficient vector.

def early_feature_importance(fitted_pipeline, feature_names):
    """Extract maximum absolute coefficient per feature."""
    coefficients = fitted_pipeline.named_steps["model"].coef_
    importance = np.abs(coefficients).max(axis=0)

    feature_importance = pd.DataFrame({
        "feature": feature_names,
        "abs_coefficient": importance,
    })
    feature_importance["omic"] = feature_importance["feature"].str.split("::", n=1).str[0]
    return feature_importance.sort_values("abs_coefficient", ascending=False)

feature_importance = early_feature_importance(early_model, X_early.columns)

TOP_N = min(100, len(feature_importance))
top_feature_usage = (
    feature_importance.head(TOP_N)["omic"]
    .value_counts()
    .rename_axis("omic")
    .reset_index(name=f"features_in_top_{TOP_N}")
)

display(top_feature_usage)
display(feature_importance.head(20))

plt.figure(figsize=(6, 4))
plt.bar(top_feature_usage["omic"], top_feature_usage[f"features_in_top_{TOP_N}"])
plt.ylabel(f"Count among top {TOP_N} coefficients")
plt.title("Which omics dominate early-integration coefficients?")
plt.xticks(rotation=30, ha="right")
plt.tight_layout()
plt.show()

	omic	features_in_top_100
0	transcriptomics	56
1	proteomics	23
2	methylation	21

	feature	abs_coefficient	omic
223955	methylation::cg26591149	0.003857	methylation
30088	proteomics::CDK1_pT14	0.003600	proteomics
30347	proteomics::PLK1	0.003514	proteomics
30119	proteomics::CYCLINB1	0.003480	proteomics
229474	methylation::cg24964883	0.003210	methylation
30136	proteomics::DNMT1	0.003144	proteomics
219956	methylation::cg16050326	0.003037	methylation
18078	transcriptomics::ENSG00000223392.1	0.002989	transcriptomics
179328	methylation::cg21250931	0.002979	methylation
30084	proteomics::cdc25C	0.002977	proteomics
30111	proteomics::Cox2	0.002966	proteomics
18200	transcriptomics::ENSG00000223949.8	0.002958	transcriptomics
30271	proteomics::MSH6	0.002941	proteomics
30038	proteomics::Aurora-A	0.002797	proteomics
29739	transcriptomics::ENSG00000287564.1	0.002743	transcriptomics
9904	transcriptomics::ENSG00000159495.8	0.002687	transcriptomics
209145	methylation::cg07980671	0.002674	methylation
5718	transcriptomics::ENSG00000127589.4	0.002652	transcriptomics
16308	transcriptomics::ENSG00000198914.5	0.002610	transcriptomics
13759	transcriptomics::ENSG00000179826.7	0.002590	transcriptomics

Interpretation (click to expand)

Proteomics has a compact feature space (464 proteins), so individual features contribute relatively more weight compared with larger omics. However, its weaker standalone performance suggests that these features capture only part of the subtype signal. The presence of proteomic features among the top-ranked predictors means that proteomics may provide complementary information rather than being uninformative.

Transcriptomics has both the strongest standalone predictive performance and the largest contribution among highly weighted features (54% of the top 100). Its dominance is therefore not explained by feature count alone; it reflects a strong subtype-associated signal distributed across many genes.

Methylation has the largest feature space (~200,000 CpG sites), but only 29% of the top 100 features are methylation-derived. This suggests that predictive methylation signal may be sparse and distributed across a small subset of sites, making it harder for simple models to identify without additional structure or feature selection.

9. Cross-omic redundancy among highly weighted features

Among the top features we just found, are different omics contributing genuinely independent information, or are some of them just correlated proxies for the same underlying signal, counted twice?

Concatenation does not model correlation between features — it has no way to know that a transcript and a protein it produces are, biologically, close to the same variable measured twice. If the top weighted features from two omics turn out to be highly correlated with each other, the model isn’t combining three independent views of the tumor; it’s partly double-counting one view. That’s a second, more specific way early integration’s flexibility can work against it (the first was Section 8’s single-omic dominance), and it further motivates keeping omics as separate, explicitly weighted models in the next section.

def top_features_per_omic(feature_importance, n_per_omic=20):
    """Return the top coefficient features within each omic."""
    selected = []
    for _, group in feature_importance.groupby("omic"):
        selected.extend(group.head(n_per_omic)["feature"].tolist())
    return selected

selected_features = top_features_per_omic(feature_importance, n_per_omic=20)
X_selected = X_early.loc[:, selected_features]

corr = X_selected.corr().abs()
rows = []

for i, feature_1 in enumerate(corr.columns):
    omic_1 = feature_1.split("::", 1)[0]
    for feature_2 in corr.columns[i + 1:]:
        omic_2 = feature_2.split("::", 1)[0]
        if omic_1 != omic_2:
            rows.append({
                "feature_1": feature_1,
                "feature_2": feature_2,
                "omic_pair": f"{omic_1} ↔ {omic_2}",
                "abs_correlation": corr.loc[feature_1, feature_2],
            })

cross_omic_corr = pd.DataFrame(rows).sort_values("abs_correlation", ascending=False)
display(cross_omic_corr.head(20))

plt.figure(figsize=(6, 4))
plt.hist(cross_omic_corr["abs_correlation"], bins=30)
plt.xlabel("Absolute correlation")
plt.ylabel("Feature-pair count")
plt.title("Cross-omic correlation among highly weighted features")
plt.tight_layout()
plt.show()

	feature_1	feature_2	omic_pair	abs_correlation
815	proteomics::CDK1_pT14	transcriptomics::ENSG00000138658.16	proteomics ↔ transcriptomics	0.439948
855	proteomics::CYCLINB1	transcriptomics::ENSG00000138658.16	proteomics ↔ transcriptomics	0.432839
835	proteomics::PLK1	transcriptomics::ENSG00000138658.16	proteomics ↔ transcriptomics	0.420436
824	proteomics::PLK1	transcriptomics::ENSG00000127589.4	proteomics ↔ transcriptomics	0.410333
1010	proteomics::ASNS	transcriptomics::ENSG00000225208.1	proteomics ↔ transcriptomics	0.403460
844	proteomics::CYCLINB1	transcriptomics::ENSG00000127589.4	proteomics ↔ transcriptomics	0.396262
1004	proteomics::ASNS	transcriptomics::ENSG00000127589.4	proteomics ↔ transcriptomics	0.393260
1075	proteomics::TRIP13	transcriptomics::ENSG00000138658.16	proteomics ↔ transcriptomics	0.392305
859	proteomics::CYCLINB1	transcriptomics::ENSG00000129484.14	proteomics ↔ transcriptomics	0.389710
804	proteomics::CDK1_pT14	transcriptomics::ENSG00000127589.4	proteomics ↔ transcriptomics	0.384738
839	proteomics::PLK1	transcriptomics::ENSG00000129484.14	proteomics ↔ transcriptomics	0.381803
970	proteomics::CYCLINE1	transcriptomics::ENSG00000225208.1	proteomics ↔ transcriptomics	0.375960
850	proteomics::CYCLINB1	transcriptomics::ENSG00000225208.1	proteomics ↔ transcriptomics	0.364153
830	proteomics::PLK1	transcriptomics::ENSG00000225208.1	proteomics ↔ transcriptomics	0.362718
819	proteomics::CDK1_pT14	transcriptomics::ENSG00000129484.14	proteomics ↔ transcriptomics	0.352746
944	proteomics::Aurora-A	transcriptomics::ENSG00000127589.4	proteomics ↔ transcriptomics	0.349084
895	proteomics::cdc25C	transcriptomics::ENSG00000138658.16	proteomics ↔ transcriptomics	0.343621
924	proteomics::MSH6	transcriptomics::ENSG00000127589.4	proteomics ↔ transcriptomics	0.340917
969	proteomics::CYCLINE1	transcriptomics::ENSG00000241359.1	proteomics ↔ transcriptomics	0.339857
1090	proteomics::FOXM1	transcriptomics::ENSG00000225208.1	proteomics ↔ transcriptomics	0.339588

Interpretation (click to expand)

The highly weighted features from early integration show moderate cross-omic correlation, concentrated mainly between proteomics and transcriptomics. This pattern is consistent with shared biological processes being captured at both the RNA and protein levels. In contrast, methylation features show less overlap with these highly weighted features, suggesting that they may provide more distinct information.

These results suggest that early integration may treat correlated features across modalities as separate sources of evidence rather than explicitly recognizing their shared biological origin. In this dataset, concatenation may therefore capture redundant signals alongside complementary ones.

→ Early integration cannot distinguish complementary signal from redundant signal automatically. To address this, we next test approaches that keep omics separate and combine their predictions through late integration.

10. Late integration by prediction averaging

Instead of merging features, what if we let each omic make its own prediction and then vote — and how much does the voting rule matter?

Late integration trains one model per omic and averages the class probabilities. Each omic gets an equal vote by default.This sidesteps the dimensionality problem of early integration and keeps each omic’s model easy to inspect on its own, but it assumes equal weighting is the right call, an assumption we will stress-test later.

def late_integrated_probabilities(probability_tables, weights=None):
    """Average same-index, same-column probability tables from omic-specific models."""
    weights = weights or {name: 1.0 for name in probability_tables}
    weight_total = sum(weights.values())

    averaged = None
    for name, proba in probability_tables.items():
        weighted_proba = proba * weights[name]
        averaged = weighted_proba if averaged is None else averaged + weighted_proba

    return averaged / weight_total

late_proba = late_integrated_probabilities(single_omic_probabilities)
late_pred = late_proba.idxmax(axis=1).rename("prediction")

results.append(evaluate_predictions(y[test_ids], late_pred, "late integration: equal weights"))
results_df = pd.DataFrame(results).sort_values("balanced_accuracy", ascending=False)
results_df

late integration: equal weights
───────────────────────────────
  Accuracy          : 0.832 
  Balanced accuracy : 0.772

	model	accuracy	balanced_accuracy
0	single omic: transcriptomics	0.824	0.798667
4	late integration: equal weights	0.832	0.772000
3	early integration: concatenation	0.752	0.702333
2	single omic: methylation	0.712	0.688333
1	single omic: proteomics	0.720	0.643333

What is your interpretation of these results?

Interpretation (click to expand)

Late integration with equal weights did not improve performance (balanced accuracy = 0.637), falling below transcriptomics alone (0.669) and even below early integration (0.647). This suggests that treating all omics as equally informative is not optimal for this dataset.

Transcriptomics provides the strongest standalone signal, while proteomics and methylation contribute weaker predictions. Equal-weight averaging assumes that each modality provides equally reliable information, so the weaker omic predictions may dilute the stronger transcriptomic signal rather than enhance it.

However, this does not mean that proteomics or methylation are uninformative. Their disagreement with transcriptomics may reflect complementary information that is not captured by a simple averaging rule. A learned integration approach could allow the model to determine how much each modality should contribute rather than imposing equal weights.

11. Late-integration weight sensitivity

How much does the (essentially arbitrary) choice of omic weights actually change the outcome — and could we have learned better weights instead?

Equal weighting is simple, but it encodes a strong assumption: each omic should contribute equally. This small test boosts one omic at a time and checks whether the held-out accuracy changes. If it does, equal weighting was an arbitrary choice, not a principled one — which is exactly what the next section addresses by learning the weights instead of guessing them.

weight_rows = []
base_weights = {name: 1.0 for name in X_omics}

# Equal weights.
proba = late_integrated_probabilities(single_omic_probabilities, weights=base_weights)
pred = proba.idxmax(axis=1)
weight_rows.append({
    "weighting": "equal weights",
    **evaluate_predictions(y[test_ids], pred, "late"),
})

# Boost each omic while keeping the same trained single-omic models and same test split.
for boosted_omic in X_omics:
    weights = base_weights.copy()
    weights[boosted_omic] = 2.0
    proba = late_integrated_probabilities(single_omic_probabilities, weights=weights)
    pred = proba.idxmax(axis=1)
    weight_rows.append({
        "weighting": f"boost {boosted_omic}",
        **evaluate_predictions(y[test_ids], pred, "late"),
    })

weight_sensitivity = (
    pd.DataFrame(weight_rows)
    .drop(columns="model")
    .sort_values("balanced_accuracy", ascending=False)
)

display(weight_sensitivity)

plt.figure(figsize=(7, 4))
plt.barh(weight_sensitivity["weighting"], weight_sensitivity["balanced_accuracy"])
plt.xlabel("Balanced accuracy on held-out test set")
plt.title("Late-integration weighting sensitivity")
plt.gca().invert_yaxis()
plt.tight_layout()
plt.show()

late
────
  Accuracy          : 0.832 
  Balanced accuracy : 0.772

late
────
  Accuracy          : 0.824 
  Balanced accuracy : 0.764

late
────
  Accuracy          : 0.832 
  Balanced accuracy : 0.707

late
────
  Accuracy          : 0.776 
  Balanced accuracy : 0.729

	weighting	accuracy	balanced_accuracy
0	equal weights	0.832	0.772000
1	boost transcriptomics	0.824	0.764000
3	boost methylation	0.776	0.729000
2	boost proteomics	0.832	0.707333

Interpretation (click to expand)

All fixed weighting schemes perform within a narrow range (balanced accuracy = 0.627–0.644), indicating that manual adjustments to modality weights produce only modest changes compared with equal weighting.

The best fixed-weight result comes from boosting proteomics (0.644), despite proteomics being the weakest standalone predictor (0.564). Conversely, boosting transcriptomics, the strongest individual modality (0.669), provides only a small improvement over equal weighting (0.641 vs 0.637), while increasing methylation contribution reduces performance (0.627).

This demonstrates that standalone predictive performance does not determine how much a modality should contribute during multimodal integration. A weaker modality may still provide complementary information when combined with stronger modalities, while a strong modality may already capture overlapping signal.

Fixed weighting is therefore a limited approach for our dataset: the optimal contribution of each modality cannot be reliably inferred from single-omic performance alone. Instead, a learned fusion strategy can adaptively determine how different modality-specific predictions should contribute to the final decision.

12. Late integration, take two — stacking (a learned combiner)

Key question: could we replace the guessed, fixed weights of equal-weight averaging with weights the data itself tells us to use?

The previous section showed that equal-weight averaging is sensitive to a choice we can’t justify from data alone. Stacking keeps the same “one model per omic” structure, but replaces the fixed averaging rule with a second linear model that learns how much each omic contributes per class. It’s still late integration in the strict sense (each omic still collapses to an independent prediction before anything cross-omic happens) but the combining rule is now learned rather than assumed.

To avoid leaking test information into the meta-learner, the per-omic probabilities used to train it are computed out-of-fold on the training set (cross-validated), not from models that have already seen the training labels. At test time we reuse the same single-omic models and probabilities from Section 6.

from sklearn.model_selection import cross_val_predict

N_FOLDS = 5

def out_of_fold_probabilities(model, X_train, y_train, class_names, n_folds=N_FOLDS):
    """Cross-validated (leakage-free) train-set probabilities for meta-learner training."""
    oof = cross_val_predict(
        clone(model), X_train, y_train, cv=n_folds, method="predict_proba"
    )
    return pd.DataFrame(oof, index=X_train.index, columns=sorted(y_train.unique()))[class_names]

meta_train_blocks = []
for name, X in X_omics.items():
    oof_proba = out_of_fold_probabilities(base_model, X.loc[train_ids], y[train_ids], class_names)
    oof_proba.columns = [f"{name}::{c}" for c in oof_proba.columns]
    meta_train_blocks.append(oof_proba)
meta_X_train = pd.concat(meta_train_blocks, axis=1)

meta_test_blocks = []
for name in X_omics:
    proba = single_omic_probabilities[name][class_names].copy()
    proba.columns = [f"{name}::{c}" for c in class_names]
    meta_test_blocks.append(proba)
meta_X_test = pd.concat(meta_test_blocks, axis=1)

meta_model = LogisticRegression(
    l1_ratio=0, C=1.0, solver="lbfgs", class_weight="balanced",
    max_iter=2000, random_state=RANDOM_STATE,
)
meta_model.fit(meta_X_train, y[train_ids])
stacked_pred = pd.Series(meta_model.predict(meta_X_test), index=test_ids, name="prediction")

results.append(evaluate_predictions(y[test_ids], stacked_pred, "late integration: stacking (meta-learner)"))
results_df = pd.DataFrame(results).sort_values("balanced_accuracy", ascending=False)
results_df

late integration: stacking (meta-learner)
─────────────────────────────────────────
  Accuracy          : 0.808 
  Balanced accuracy : 0.801

	model	accuracy	balanced_accuracy
5	late integration: stacking (meta-learner)	0.808	0.801333
0	single omic: transcriptomics	0.824	0.798667
4	late integration: equal weights	0.832	0.772000
3	early integration: concatenation	0.752	0.702333
2	single omic: methylation	0.712	0.688333
1	single omic: proteomics	0.720	0.643333

# Unlike late integration's fixed weights, stacking's learned coefficients show
# how much the meta-learner trusts each omic, per predicted subtype.
meta_coefs = pd.DataFrame(
    meta_model.coef_, index=meta_model.classes_, columns=meta_X_train.columns
)
meta_coefs.columns = pd.MultiIndex.from_tuples(
    [tuple(col.split("::")) for col in meta_coefs.columns], names=["omic", "class"]
)

#magnitude of reliance
meta_coefs.abs().T.groupby(level="omic").mean().T.round(2)

omic	methylation	proteomics	transcriptomics
Basal	0.36	1.02	0.75
Her2	0.50	0.62	0.92
LumA	0.49	0.49	0.72
LumB	0.73	0.57	0.76
Normal	0.52	0.44	1.13

Interpretation (click to expand)

The stacking model does not assign one fixed contribution to each omic. Instead, the meta-model learns subtype-specific contributions for each modality prediction, allowing different omics to influence different classification decisions.

For example, transcriptomics receives larger coefficients for Her2 and Normal predictions compared with methylation, while its contribution is closer to proteomics for Basal and LumA. This type of class-specific flexibility cannot be captured by equal-weight late integration, which assumes that all modalities contribute equally across all predictions.

These coefficients illustrate why learned fusion can be advantageous: the optimal contribution of each modality depends on the prediction context rather than being determined solely by standalone performance. However, these coefficients describe prediction combination rather than direct biological importance of each omic.

13. Compare all linear methods

After all that, did integration actually beat the best single omic? And which integration strategy “won” ?

A genuine improvement here justifies the added complexity of early, late, or stacked integration. A tie (or a loss) can be just as informative; it tells you that the extra modeling complexity is not justified for this dataset.

summary = pd.DataFrame(results).sort_values("balanced_accuracy", ascending=True)
display(summary.sort_values("balanced_accuracy", ascending=False))

plt.figure(figsize=(8, 4.5))
plt.barh(summary["model"], summary["balanced_accuracy"])
plt.xlabel("Balanced accuracy on held-out test set")
plt.title("Single-omic vs integrated linear models")
plt.tight_layout()
plt.show()

	model	accuracy	balanced_accuracy
5	late integration: stacking (meta-learner)	0.808	0.801333
0	single omic: transcriptomics	0.824	0.798667
4	late integration: equal weights	0.832	0.772000
3	early integration: concatenation	0.752	0.702333
2	single omic: methylation	0.712	0.688333
1	single omic: proteomics	0.720	0.643333

Interpretation (click to expand)

Naive integration strategies did not improve over the strongest single-omic model in this benchmark. Early integration through concatenation (0.647 balanced accuracy) and late integration through equal-weight averaging (0.637) both performed below transcriptomics alone (0.669), suggesting that simply combining modalities does not guarantee improved prediction.

In contrast, stacking achieved the highest performance (0.699 balanced accuracy), exceeding transcriptomics alone by learning how to combine modality-specific predictions rather than applying a fixed fusion rule. The learned combination allows different omics to contribute differently depending on the classification context, capturing complementary information that was not fully exploited by simpler integration strategies.

This result demonstrates that the challenge in multi-omics integration is not simply adding more data, but learning how different modalities should contribute to the final prediction.

14. Confusion matrix and classification report

Use this section to discuss which subtypes are easier or harder for the best simple model.

best_model_name = results_df.iloc[0]["model"]
print(f"Best model by balanced accuracy: {best_model_name}")

if best_model_name == "early integration: concatenation":
    best_pred = early_pred
elif best_model_name == "late integration: equal weights":
    best_pred = late_pred
elif best_model_name == "late integration: stacking (meta-learner)":
    best_pred = stacked_pred
else:
    best_omic = best_model_name.split(": ", 1)[1]
    best_pred = single_omic_models[best_omic].predict(X_omics[best_omic].loc[test_ids])
    best_pred = pd.Series(best_pred, index=test_ids)
    
plot_confusion_matrix(y[test_ids], best_pred)

Best model by balanced accuracy: late integration: stacking (meta-learner)

Interpretation (click to expand)

The stacked model separates Basal and Her2 very well, with near-perfect recall for both classes. This suggests that these subtypes contain strong predictive signals that are consistently captured across the integrated modality predictions.

The main challenge is distinguishing LumA from LumB. These subtypes show substantial overlap, with many LumB samples classified as LumA, consistent with the idea that these molecular subtypes represent a biological continuum rather than completely distinct categories.

The Normal class shows high precision but low recall: when the model predicts Normal it is usually correct, but it misses many true Normal samples. This likely reflects both the limited number of Normal samples available (only 7 test patients) and the difficulty of learning a minority class with few examples.

Overall, this demonstrates that integration strategy can improve how information is combined across omics, but it cannot compensate for limited sample size, class imbalance, or intrinsically overlapping biological categories.

Summary

Strategy	Data Flow	Integration Point	Output	Key Characteristic
Early Integration	RNA-seq + Methylation + Proteomics → Combine	Before training	One feature matrix	All modalities are merged into a single input representation.
Late Integration	RNA-seq → Model Methylation → Model Proteomics → Model → Average predictions	After prediction	Combined prediction	Each modality is modeled independently before combining predictions.
Stacking	RNA-seq → Model Methylation → Model Proteomics → Model → Meta-model	Meta-learning	Final prediction	Learns the optimal way to combine modality-specific model outputs.

Takeaways

1. Different omics capture different aspects of biological variation. PCA provides a visual entry point for comparing how modalities represent sample structure and where shared or distinct patterns may exist.

2. Single-omic models are essential baselines. Integration should be evaluated against individual modalities to determine whether combining data provides additional value beyond single sources.

3. Early integration is simple but assumes features can be combined directly. Concatenation creates a unified feature space, allowing models to use information across omics simultaneously. However, it does not explicitly model relationships between modalities, so correlated features may be treated as independent evidence.

4. Predictive features are not necessarily independent features. Correlated features across omics can represent shared biological processes, meaning that high feature importance does not automatically indicate unique information.

5. Late integration is interpretable but depends on the fusion rule. Combining modality-specific predictions keeps models separate and easier to interpret, but fixed weighting schemes may not capture differences in modality reliability.

6. Stacking provides a more adaptive form of late integration. Instead of fixing how modalities contribute, a meta-model learns how to combine modality-specific predictions from the data.

7. Classical integration methods do not explicitly model cross-omic feature relationships. This limitation motivates approaches that learn shared structure across modalities, covered next in Part 2.

8. The goal is not to find the most complex integration method, but the method whose assumptions match the data and biological question.

Bridge to the next part

None of the classical integration methods we just learned explicitly model relationships between features across omics. Feature concatenation places all features into a shared space but does not account for redundancy or structured relationships between modalities (which is why proteomics/transcriptomics redundancy in Section 9 was not apparent until we examined feature relationships directly). Prediction averaging and stacking address this limitation indirectly; they combine modality-specific predictions without examining the feature-level relationships that generated them.

A method that models cross-omic relationships directly could separate shared signal (such as the genes we found correlated across proteomics and transcriptomics) from genuinely omic-specific variation. This avoids treating redundant features as independent evidence during early integration or relying only on prediction-level combinations during stacking. That is precisely what the next session introduces: methods that learn shared latent factors directly from cross-omic covariance, making the redundancy we identified manually an explicit, modeled part of the representation rather than a blind spot.

These approaches are not necessarily “better” than classical integration methods; they make different assumptions about where useful information exists in the data. Their advantage comes from being able to model shared and modality-specific structure explicitly when those relationships are important for the biological question.

This moves us towards generating an integrated multi-omic profile for each patient, where shared and modality-specific biological signals can be represented explicitly.