Part 1 — MOFA tools for an LLM agent (LangChain `@tool`)#

In this notebook, we’re going to learn about how to make LLMs use tools. Tools are simply functions, given to a LLM agent, which then decides which to call to answer a given task, instead of relying only on its internal knowledge.

We’ll specifically be using as tools the functions defined in src/mofa_tools.py. These tools will help Claude answer questions about the multi-omics MOFA model we fitted in Part 0: A fitted MOFA model of TCGA breast-cancer subtypes, with 603 patients × 3 omics views (transcriptomics, proteomics, methylation), PAM50 subtype per patient.

For Claude to know when to use each function, they all contain a short description (a docstring) of what it does, and Claude decides on its own which ones to use, and in what order.

Key design choice: the model is not fit by the agent. Within these functions, we have decided not to include the function to fit the MOFA model (fit_mofa). The agent’s tools operate on the already fitted model. Fitting is expensive and non-interactive, so it is deliberately not exposed as a tool. Tools are the cheap, interpretable analysis steps.

Learning objectives#

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

Explain what a tool is, and why handing Claude a function is not the same as trusting it with what is inside that function.
Write a plain Python function that returns tool-friendly evidence, and wrap it as a tool Claude can call using LangChain’s @tool decorator.
Write a system prompt that gives Claude the grounding it needs to answer questions safely.
Read an agent’s tool-call trace and judge whether it combined evidence from more than one tool correctly.

0. Environment & API key#

Load the key, and point the notebook at the project root so it can find src/mofa_tools.py.

from pathlib import Path
import os, sys, json

from dotenv import load_dotenv

SESSION_DIR = Path(os.environ.get(
    "ECCB_SESSION3_DIR",
    Path.home() / "ECCB2026_TEST/sessions/session-3-agentic-llm-workflows"))
assert SESSION_DIR.exists(), f"Session directory not found: {SESSION_DIR}"
DOTENV_PATH = "/.env"
DATA_DIR = Path("/data") if Path("/data").exists() else SESSION_DIR / "data"  # shared mount on workshop servers

PROJECT_ROOT = SESSION_DIR            # anchors imports (src/...) and paths below
sys.path.insert(0, str(PROJECT_ROOT))
load_dotenv(DOTENV_PATH)

if not os.environ.get("ANTHROPIC_API_KEY"):
    raise RuntimeError("ANTHROPIC_API_KEY not found. Add it to a .env in the project root.")
print("ANTHROPIC_API_KEY loaded:", bool(os.environ.get("ANTHROPIC_API_KEY")))

ANTHROPIC_API_KEY loaded: True

1. Load the data + the cached MOFA model (one-time setup)#

This mirrors the pipeline in src/mofa_tools.py: load aligned omics, make the shared train/test split, then load the pre-fitted MOFA model from .hdf5 (fit once via the project’s fit step, never re-fit here). The agent’s tools close over the objects built in this cell.

import numpy as np
import pandas as pd
import mofax as mfx

from src.mofa_tools import (
    RANDOM_STATE, TEST_SIZE, N_TOP_VARIABLE_HIGH_DIM_FEATURES,
    HIGH_DIMENSIONAL_VIEWS, MAX_FACTORS, MIN_TOTAL_R2,
    load_omics_data, make_train_test_split, select_active_factors,
    project_test_patients_to_mofa_factors, eta_squared_by_factor,
    fit_factor_classifier,
)

CACHE = PROJECT_ROOT / "outputs" / (
    f"trained_mofaplus_train_var{N_TOP_VARIABLE_HIGH_DIM_FEATURES}"
    f"_max{MAX_FACTORS}_ard_model.hdf5")
assert CACHE.exists(), f"Cached MOFA model not found: {CACHE}. Run the one-time fit first."

X_omics, y = load_omics_data(DATA_DIR)
X_train, X_test, y_train, y_test, train_ids, test_ids = make_train_test_split(
    X_omics, y, TEST_SIZE, RANDOM_STATE, HIGH_DIMENSIONAL_VIEWS,
    N_TOP_VARIABLE_HIGH_DIM_FEATURES)

view_names = list(X_train.keys())
train_ids_str, test_ids_str = train_ids.astype(str), test_ids.astype(str)

model = mfx.mofa_model(str(CACHE))                       # LOAD, never fit
train_factors = model.get_factors(df=True)
train_factors.index = train_factors.index.astype(str)
test_factors = project_test_patients_to_mofa_factors(
    model, X_train, X_test, train_factors, view_names)
factors_df = pd.concat([train_factors, test_factors], axis=0)
factors_df.columns = factors_df.columns.astype(str)
factors_df = factors_df.reindex(y.index.astype(str))

active_cols, r2_summary = select_active_factors(model, MIN_TOTAL_R2, MAX_FACTORS)
r2_all = model.get_r2().rename(columns={"Factor": "factor", "View": "view",
                                        "Group": "group", "R2": "r2"})
print(f"Loaded cached MOFA model | views={view_names}")
print(f"train={len(train_ids)}  test={len(test_ids)}  active factors={active_cols}")

Loaded cached MOFA model | views=['transcriptomics', 'proteomics', 'methylation']
train=375  test=125  active factors=['Factor1', 'Factor2', 'Factor3', 'Factor4', 'Factor5', 'Factor6', 'Factor7', 'Factor8', 'Factor9', 'Factor10']

2. The tool functions (plain Python, on the cached model)#

Each function returns a small JSON-serializable result, the kind of compact, grounded evidence an LLM can reason over. We demonstrate each directly here before handing them to the model.

#################################################
# Write a funciton & its docstring that summarises 
# the patients information and features per omics view, plus the subtype class counts.
def data_summary():
    """Patients and features per omics view, plus the subtype class counts."""
    return {"views": {v: {"patients": int(X_omics[v].shape[0]),
                          "features": int(X_omics[v].shape[1])} for v in view_names},
            "subtype_counts": y.value_counts().to_dict()}
#################################################

def split_summary():
    """Train/test patient counts and the number of features kept per view after
    variable-feature selection."""
    return {"n_train": int(len(train_ids)), "n_test": int(len(test_ids)),
            "features_per_view": {v: int(X_train[v].shape[1]) for v in view_names}}

def active_factors():
    """Which MOFA factors are 'active' (total variance explained >= threshold),
    with each factor's total R2 summed across views."""
    return {"active_factors": active_cols,
            "total_r2": r2_summary.set_index("factor")["total_r2"].round(3).to_dict()}

def factor_view_r2(factor):
    """Variance explained (R2, %) by one factor in each omics view, and the view
    that factor explains most."""
    sub = r2_all[r2_all["factor"] == factor][["view", "r2"]]
    d = {k: round(float(v), 3) for k, v in zip(sub["view"], sub["r2"])}
    return {"factor": factor, "r2_by_view": d,
            "top_view": max(d, key=d.get) if d else None}

def factor_subtype_association():
    """Rank active factors by eta-squared: how much of each factor's variance is
    explained by PAM50 subtype (training patients). Higher = more subtype-linked."""
    assoc = eta_squared_by_factor(factors_df.loc[train_ids_str, active_cols], y_train)
    return assoc.assign(eta_squared=assoc["eta_squared"].round(3)).to_dict("records")

def top_features_for_factor(factor, view="transcriptomics", n=5):
    """Top positive- and negative-weighted features of `view` for `factor` — the
    features that most strongly define that latent axis."""
    w = model.get_weights(views=view, df=True); w.columns = w.columns.astype(str)
    s = w[factor].sort_values()
    return {"factor": factor, "view": view,
            "top_negative": {k: round(float(v), 3) for k, v in s.head(n).items()},
            "top_positive": {k: round(float(v), 3) for k, v in s.tail(n).items()}}

def classify_subtype_from_factors():
    """Fit logistic regression on MOFA factors (train) and evaluate on held-out
    patients: accuracy / balanced-accuracy / macro-F1 and the most-confused pair."""
    _, pred, metrics = fit_factor_classifier(
        factors_df[active_cols], y, train_ids_str, test_ids_str, "MOFA factors + LR")
    labels = sorted(set(y_test.values) | set(pred))
    cm = pd.crosstab(pd.Series(y_test.values, name="true"),
                     pd.Series(pred, name="pred")).reindex(
        index=labels, columns=labels, fill_value=0)
    off = cm.to_numpy(copy=True); np.fill_diagonal(off, 0)
    r, c = np.unravel_index(off.argmax(), off.shape)
    return {"metrics": {k: round(v, 3) for k, v in metrics.items() if k != "model"},
            "most_confused": {"true": labels[r], "predicted": labels[c], "count": int(off.max())}}

def train_vs_test_subtype_association():
    """Compare factor<->subtype eta-squared on TRAIN vs projected TEST patients —
    does the factor space still separate subtypes out of sample?"""
    tr = eta_squared_by_factor(factors_df.loc[train_ids_str, active_cols], y_train)
    te = eta_squared_by_factor(factors_df.loc[test_ids_str, active_cols], y_test)
    m = tr.merge(te, on="factor", suffixes=("_train", "_test"))
    m["eta_squared_train"] = m["eta_squared_train"].round(3)
    m["eta_squared_test"] = m["eta_squared_test"].round(3)
    return m.to_dict("records")

Quick direct demo — the evidence the agent will be working from:

print("data_summary          :", json.dumps(data_summary()["subtype_counts"]))
print("split_summary          :", json.dumps(split_summary()))
print("active_factors         :", active_factors()["active_factors"])
print("factor_view_r2(Factor1):", json.dumps(factor_view_r2("Factor1")))
print("subtype assoc (top 3)  :", factor_subtype_association()[:3])
print("classify               :", json.dumps(classify_subtype_from_factors()))

data_summary          : {"LumA": 237, "LumB": 100, "Basal": 97, "Her2": 41, "Normal": 25}
split_summary          : {"n_train": 375, "n_test": 125, "features_per_view": {"transcriptomics": 2000, "proteomics": 464, "methylation": 2000}}
active_factors         : ['Factor1', 'Factor2', 'Factor3', 'Factor4', 'Factor5', 'Factor6', 'Factor7', 'Factor8', 'Factor9', 'Factor10']
factor_view_r2(Factor1): {"factor": "Factor1", "r2_by_view": {"transcriptomics": 1.426, "proteomics": 0.811, "methylation": 32.407}, "top_view": "methylation"}
subtype assoc (top 3)  : [{'factor': 'Factor2', 'eta_squared': 0.824}, {'factor': 'Factor6', 'eta_squared': 0.373}, {'factor': 'Factor1', 'eta_squared': 0.281}]
classify               : {"metrics": {"accuracy": 0.776, "balanced_accuracy": 0.814, "macro_f1": 0.751}, "most_confused": {"true": "LumA", "predicted": "LumB", "count": 10}}

3. Wrap the functions as LangChain tools#

LangChain’s @tool marker turns a plain Python function into something Claude can call directly. The docstring becomes the description the model reads to decide when a tool is useful — Claude never sees the function’s actual code, only its name, arguments, and docstring. Note we do not wrap fit_mofa, build_mofa_matrix_input, or evaluate_predictions, the first is expensive setup, the other two are internal plumbing. Exposing them would only give the model ways to misfire.

from langchain_core.tools import tool

@tool
def data_summary_tool() -> dict:
    """Return patients and features per omics view, plus PAM50 subtype class counts."""
    return data_summary()

@tool
def split_summary_tool() -> dict:
    """Return train/test patient counts and features kept per view after selection."""
    return split_summary()

@tool
def active_factors_tool() -> dict:
    """Return which MOFA factors are active and each factor's total R2 across views."""
    return active_factors()

@tool
def factor_view_r2_tool(factor: str) -> dict:
    """Variance explained (R2) by one factor (e.g. 'Factor1') in each omics view,
    and which view it explains most."""
    #################################################
    return factor_view_r2(factor)
    #################################################
    
@tool
def factor_subtype_association_tool() -> list:
    """Rank active factors by eta-squared association with PAM50 subtype (train)."""
    return factor_subtype_association()

@tool
def top_features_for_factor_tool(factor: str, view: str = "transcriptomics", n: int = 5) -> dict:
    """Top positive/negative weighted features of a view for a factor (its drivers)."""
    return top_features_for_factor(factor, view=view, n=n)

@tool
def classify_subtype_from_factors_tool() -> dict:
    """Predict subtype from MOFA factors; return held-out metrics + most-confused pair."""
    return classify_subtype_from_factors()

@tool
def train_vs_test_subtype_association_tool() -> list:
    """Compare factor<->subtype eta-squared on train vs projected test patients."""
    return train_vs_test_subtype_association()

tools = [data_summary_tool, split_summary_tool, active_factors_tool, factor_view_r2_tool,
         factor_subtype_association_tool, top_features_for_factor_tool,
         classify_subtype_from_factors_tool, train_vs_test_subtype_association_tool]
tools_by_name = {t.name: t for t in tools}
list(tools_by_name)

['data_summary_tool',
 'split_summary_tool',
 'active_factors_tool',
 'factor_view_r2_tool',
 'factor_subtype_association_tool',
 'top_features_for_factor_tool',
 'classify_subtype_from_factors_tool',
 'train_vs_test_subtype_association_tool']

4. Bind the tools to Claude#

bind_tools tells the model the tool schemas: what these tools are called, what arguments they take, and what their docstrings say. Claude decides at question time whether to use any of them. Creating the client and binding it to the tools does not call the API by itself; no tokens are spent until Claude actually answers a question in Section 6.

from langchain_anthropic import ChatAnthropic

MODEL = "claude-haiku-4-5"   # swap for any current Claude model id
llm = ChatAnthropic(model=MODEL, temperature=0)
llm_with_tools = llm.bind_tools(tools)
print("Bound", len(tools), "tools to", MODEL)

Bound 8 tools to claude-haiku-4-5

5. The agent loop#

A tool-calling agent is just a loop: send the conversation so far, run whatever tools the model asks for, add the results back into the conversation, and repeat until the model answers without requesting any more tools.

from langchain_core.messages import AIMessage, HumanMessage, SystemMessage, ToolMessage

#################################################
# Write the prompt to define the persona, the activity that the agent should perform. 
# Additional clue: List the factors and the subtypes in the prompt.
SYSTEM = ("You are a computational-biology assistant analysing a fitted MOFA model "
          "of TCGA breast-cancer multi-omics data. Use the tools to gather evidence; "
          "ground every quantitative claim in tool results and name which tool you used. "
          "Factors are named 'Factor1'..'Factor10'; subtypes are PAM50 (LumA, LumB, "
          "Basal, Her2, Normal).")
#################################################

def run_agent(question: str, max_steps: int = 6, verbose: bool = True) -> AIMessage:
    messages = [SystemMessage(content=SYSTEM), HumanMessage(content=question)]
    for step in range(max_steps):

        # Ask the model what to do next
        model_response = llm_with_tools.invoke(messages)
        messages.append(model_response)

        # If no tool was requested, this is the final answer
        if not model_response.tool_calls:
            return model_response
        
        # Otherwise, run each requested tool
        for call in model_response.tool_calls:
            if verbose:
                print(f"[step {step}] -> {call['name']}({call['args']})")
            result = tools_by_name[call["name"]].invoke(call["args"])

            # Add the tool result to the conversation
            messages.append(ToolMessage(content=json.dumps(result, default=str),
                                        tool_call_id=call["id"]))
    raise RuntimeError(f"Agent did not finish within {max_steps} steps.")

6. Run the agent (these cells spend API tokens)#

Each run_agent(...) call below issues real Anthropic API requests. They go from simple, single-tool questions to ones that need Claude to combine evidence from more than one tool. Run the ones you want.

Simple, single-tool questions#

# Q1: data shape
print(run_agent("How many patients and features are in each omics view, "
                "and how many patients per PAM50 subtype?").content)

[step 0] -> data_summary_tool({})
Based on the **data_summary_tool** results, here's the breakdown:

## Patients and Features per Omics View

| View | Patients | Features |
|------|----------|----------|
| **Transcriptomics** | 500 | 29,995 |
| **Proteomics** | 500 | 464 |
| **Methylation** | 500 | 200,000 |

## Patients per PAM50 Subtype

| Subtype | Count |
|---------|-------|
| **LumA** | 237 |
| **LumB** | 100 |
| **Basal** | 97 |
| **Her2** | 41 |
| **Normal** | 25 |
| **Total** | 500 |

The dataset is dominated by Luminal A (LumA) subtype (47.4% of patients), followed by Luminal B and Basal subtypes (~20% each), with Her2 and Normal subtypes being less common (8.2% and 5% respectively).

# Q2: the split
#################################################
# Write the text to ask the agent to split the cohort between test and train. 
# Additional clue: detail the information that you are asking for, eg. how many patients each?
print(run_agent(
"Split the cohort into train/test. How many patients in each, "
"and how many features survive per view after selection?"
#################################################
).content)

[step 0] -> split_summary_tool({})
Based on the **split_summary_tool** results:

**Patient counts:**
- **Training set:** 375 patients
- **Test set:** 125 patients
- **Total:** 500 patients (75/25 split)

**Features retained per view after selection:**
- **Transcriptomics:** 2,000 features
- **Proteomics:** 464 features
- **Methylation:** 2,000 features

The feature selection has been most stringent for proteomics (464 features), while transcriptomics and methylation retained 2,000 features each.

# Q3: active factors
print(run_agent("Which MOFA factors are active, and how much total variance "
                "does each explain?").content)

[step 0] -> active_factors_tool({})
Based on the **active_factors_tool**, all 10 factors are active in the MOFA model. Here's the total variance (R²) explained by each:

| Factor | Total R² (%) |
|--------|-------------|
| **Factor1** | 34.64% |
| **Factor2** | 29.12% |
| **Factor3** | 17.22% |
| **Factor4** | 10.65% |
| **Factor5** | 8.87% |
| **Factor6** | 6.50% |
| **Factor7** | 5.93% |
| **Factor8** | 4.95% |
| **Factor9** | 4.63% |
| **Factor10** | 3.98% |

**Key observations:**
- **Factor1 and Factor2 dominate**, together explaining ~63.8% of total variance
- The first three factors account for ~81% of variance
- Variance drops off substantially after Factor3, with each subsequent factor explaining <11%
- All 10 factors combined explain 126.6% (note: this is cumulative R² across views, so values can exceed 100%)

# Q4: factor vs view
print(run_agent("Which omics view does Factor1 explain the most?").content)

[step 0] -> factor_view_r2_tool({'factor': 'Factor1'})
Based on the **factor_view_r2_tool** results, **Factor1 explains the most variance in the methylation view**, with an R² of **32.41%**. 

For comparison, Factor1 explains:
- **1.43%** of variance in transcriptomics
- **0.81%** of variance in proteomics

This indicates that Factor1 is primarily a methylation-driven factor.

Multi-step reasoning#

# Q5: FLAGSHIP = strongest subtype factor + its drivers
print(run_agent("Which MOFA factor is most associated with breast-cancer subtype, "
                "and which transcriptomic features most strongly drive it?").content)

[step 0] -> factor_subtype_association_tool({})
[step 1] -> top_features_for_factor_tool({'factor': 'Factor2', 'view': 'transcriptomics', 'n': 10})
## Summary

**Factor2** is the MOFA factor most strongly associated with breast-cancer PAM50 subtype (η² = 0.824, from `factor_subtype_association_tool`).

### Top Transcriptomic Drivers of Factor2:

**Positive drivers** (highest weights):
1. **ENSG00000160182.3** (weight: 0.589)
2. **ENSG00000173467.9** (weight: 0.568)
3. **ENSG00000160180.15** (weight: 0.532)
4. **ENSG00000106541.12** (weight: 0.493)
5. **ENSG00000235687.9** (weight: 0.492)

**Negative drivers** (lowest weights):
1. **ENSG00000166535.20** (weight: -0.365)
2. **ENSG00000185686.18** (weight: -0.36)
3. **ENSG00000164434.12** (weight: -0.352)
4. **ENSG00000186832.9** (weight: -0.349)
5. **ENSG00000102243.13** (weight: -0.344)

These transcriptomic features represent the key molecular drivers that distinguish breast-cancer subtypes in this MOFA model. The strong association (η² = 0.824) indicates that Factor2 captures substantial subtype-specific transcriptomic variation.

# Q6: predict + failure mode
print(run_agent("Can we predict subtype from the MOFA factors alone? Report the "
                "held-out performance and say which two subtypes are most confused.").content)

[step 0] -> classify_subtype_from_factors_tool({})
Yes, we can predict PAM50 subtype from MOFA factors alone with reasonable performance on held-out test data:

**Held-out Performance:**
- **Accuracy:** 77.6%
- **Balanced Accuracy:** 81.4%
- **Macro F1:** 0.751

**Most Confused Subtype Pair:**
- **LumA** (true) is most frequently misclassified as **LumB** (10 instances)

This makes biological sense: Luminal A and Luminal B are both luminal subtypes and share many molecular features, so they are naturally more difficult to distinguish than, say, Basal or Her2 subtypes. The model achieves good overall discrimination, suggesting that the MOFA factors capture subtype-relevant variation across the multi-omics data.

# Q7: variance vs subtype separation (synthesises two tools)
#################################################
# Write the text to ask the agent to compare the R2 ranking with the subtype association ranking and explain 
# any mismatch, in order to understand the relationship between variance and subtype separation
print(run_agent("Is the factor that explains the most variance also the one that "
                "best separates subtypes? Compare the R2 ranking with the subtype "
                "association ranking and explain any mismatch.").content)
#################################################

[step 0] -> active_factors_tool({})
[step 0] -> factor_subtype_association_tool({})
Excellent! Now I can provide a detailed comparison. Here's what the data shows:

## R2 Ranking (Variance Explained - Total Across All Views)
1. **Factor1**: 34.64%
2. **Factor2**: 29.12%
3. **Factor3**: 17.22%
4. **Factor4**: 10.65%
5. **Factor5**: 8.87%
6. **Factor6**: 6.50%
7. **Factor7**: 5.93%
8. **Factor8**: 4.95%
9. **Factor9**: 4.63%
10. **Factor10**: 3.98%

## Subtype Association Ranking (η² with PAM50 Subtype)
1. **Factor2**: 0.824
2. **Factor6**: 0.373
3. **Factor1**: 0.281
4. **Factor4**: 0.181
5. **Factor8**: 0.154
6. **Factor9**: 0.086
7. **Factor10**: 0.080
8. **Factor3**: 0.049
9. **Factor5**: 0.046
10. **Factor7**: 0.018

## Key Findings & Mismatch Explanation

**The factor that explains the most variance (Factor1) is NOT the one that best separates subtypes (Factor2).**

### The Mismatch:
- **Factor1** dominates in total variance (34.64% R²) but ranks **3rd** in subtype association (η² = 0.281)
- **Factor2** is 2nd in variance (29.12% R²) but ranks **1st** in subtype association (η² = 0.824)

### Why This Mismatch Occurs:

1. **Factor1 captures general biological variation**: With the highest R² across all omics views, Factor1 likely captures broad, subtype-independent biological signals (e.g., general tumor vs. normal tissue differences, or general proliferation signals that span multiple subtypes).

2. **Factor2 is subtype-specific**: Despite explaining less total variance, Factor2 has a dramatically higher association with PAM50 subtypes (η² = 0.824 vs. 0.281). This suggests Factor2 captures variation that is highly discriminative between the five PAM50 classes, even if it doesn't explain as much overall variance.

3. **Factor6 is another interesting case**: It ranks 6th in variance (6.50% R²) but 2nd in subtype association (η² = 0.373), further illustrating that subtype-discriminative power is independent of total variance explained.

### Biological Interpretation:
This is actually a **healthy pattern** in multi-omics analysis. The largest factors often capture ubiquitous biological signals (like proliferation or immune infiltration), while smaller, more focused factors capture the specific molecular signatures that define disease subtypes. Factor2 appears to be precisely such a subtype-defining factor.

# Q8: out-of-sample check
print(run_agent("Do the MOFA factors still separate subtypes on the held-out test "
                "patients? Compare the train and test subtype associations.").content)

[step 0] -> train_vs_test_subtype_association_tool({})
Excellent news: **Yes, the MOFA factors maintain or even strengthen their ability to separate subtypes on held-out test patients.** Here's the detailed comparison:

## Train vs Test Subtype Association (η² = effect size)

| Factor | Train η² | Test η² | Change |
|--------|----------|---------|--------|
| **Factor2** | 0.824 | 0.872 | ↑ +0.048 |
| **Factor6** | 0.373 | 0.365 | ↓ -0.008 |
| **Factor1** | 0.281 | 0.399 | ↑ +0.118 |
| **Factor4** | 0.181 | 0.277 | ↑ +0.096 |
| **Factor8** | 0.154 | 0.299 | ↑ +0.145 |
| **Factor9** | 0.086 | 0.115 | ↑ +0.029 |
| **Factor10** | 0.080 | 0.103 | ↑ +0.023 |
| **Factor3** | 0.049 | 0.181 | ↑ +0.132 |
| **Factor5** | 0.046 | 0.097 | ↑ +0.051 |
| **Factor7** | 0.018 | 0.091 | ↑ +0.073 |

## Key Findings:

1. **Strong generalization**: 9 out of 10 factors show *improved* or stable subtype associations on test data
2. **Factor2 dominates**: Maintains the strongest association (η² = 0.872), explaining ~87% of subtype variance
3. **Surprising improvements**: Several factors (Factor8, Factor3, Factor1) show notably *stronger* associations on test data, suggesting the model captures robust biological signals
4. **No overfitting**: The lack of degradation indicates the factors are not overfitted to training data

This suggests the MOFA model has learned generalizable patterns that robustly separate PAM50 subtypes across independent patient cohorts.

Reflection#

Which tool calls did the model choose for each question, and in what order?
For Q5/Q7, did it combine evidence from multiple tools, or stop early?
What would be unsafe for the model to state without calling these tools?
We intentionally did not expose fit_mofa. What goes wrong if an agent can trigger an expensive, non-deterministic fit mid-conversation?
The MCP track (02_MCP_v1.ipynb) and skills track (03_agent_skills_v1.ipynb) reuse this exact mofa_tools.py backend, compare how the tool layer is exposed.

Part 1 — MOFA tools for an LLM agent (LangChain @tool)#