Causal Discovery on Dependent Mixed Data with Applications to Gene Regulatory Network Inference

1. Problemstellung

Das Paper adressiert zwei fundamentale Herausforderungen bei der kausalen Entdeckung (Causal Discovery) aus Beobachtungsdaten:

1. Verletzung der i.i.d.-Annahme: Die meisten bestehenden Methoden zur Inferenz von Directed Acyclic Graphs (DAGs) gehen davon aus, dass die Beobachtungseinheiten (Samples) unabhängig und identisch verteilt (i.i.d.) sind. In realen Anwendungen wie sozialen Netzwerken oder der Genomik (z. B. Einzelzell-RNA-Sequenzierung, scRNA-seq) sind Proben jedoch oft voneinander abhängig (z. B. durch gemeinsame Abstammungslinien, Batch-Effekte oder räumliche Organisation). Diese Abhängigkeiten führen zu Korrelationen zwischen den Einheiten, die die Schätzung kausaler Strukturen verzerren.
2. Gemischte Datentypen (Mixed Data): Viele Datensätze enthalten eine Mischung aus kontinuierlichen Variablen (z. B. Genexpressionsniveaus) und diskreten Variablen (z. B. binäre Zustände oder kategoriale Stufen). Die Modellierung von Abhängigkeiten zwischen Proben ist bei diskreten Daten besonders schwierig, da Standardmethoden für kontinuierliche Daten nicht direkt anwendbar sind.

Ziel ist es, einen DAG zu lernen, der die kausalen Beziehungen zwischen Variablen korrekt erfasst, obwohl die Daten sowohl gemischte Typen aufweisen als auch Abhängigkeiten zwischen den Beobachtungseinheiten bestehen.

2. Methodik

Die Autoren schlagen ein De-Korrelations-Framework vor, das auf einem strukturellen Gleichungsmodell (SEM) mit latenten Variablen basiert. Der Ansatz gliedert sich in folgende Hauptkomponenten:

A. Das Modell für abhängige gemischte Daten

Das Modell erweitert bestehende SEMs durch die Einführung latenter kontinuierlicher Variablen $Z$:

• Latente Struktur: Für jede beobachtete Variable $X_j$ wird eine latente Variable $Z_j$ angenommen.
• Generierungsprozess:
  • Für kontinuierliche Variablen ($j \in \mathcal{C}$): $X_j = Z_j$.
  • Für diskrete Variablen ($j \in \mathcal{D}$): $X_j$ entsteht durch Quantisierung (Thresholding) von $Z_j$ mittels Schwellenwerten $T_j$.
  • Die latente Variable folgt einem SEM: $Z_j = \sum_{k \in PA_j} \beta_{kj} g(Z_k) + \varepsilon_j$, wobei $g(\cdot)$ die Identitätsfunktion für kontinuierliche Eltern und die Quantisierungsfunktion für diskrete Eltern ist.
• Abhängigkeitsmodellierung: Der Fehlerterm $\varepsilon_j$ folgt einer multivariaten Normalverteilung $N_n(0, \Sigma)$. Die Kovarianzmatrix $\Sigma$ modelliert die Abhängigkeiten zwischen den Einheiten (Samples), nicht zwischen den Variablen. Dies induziert Korrelationen zwischen den Zeilen des Datenmatrizen.

B. De-Korrelations-Strategie (De-correlation)

Das Kernkonzept besteht darin, die Abhängigkeiten zwischen den Einheiten zu entfernen, um Standard-Algorithmen für kausale Entdeckung anwenden zu können.

• Präzisionsmatrix: Sei $\Theta = \Sigma^{-1}$ die Präzisionsmatrix und $L$ deren Cholesky-Zerlegung ($\Theta = LL^\top$).
• Transformation: Durch Anwendung von $L^\top$ auf die Fehlerterme (und damit auf die latenten Daten) werden die Einheiten de-korreliert: $\tilde{\varepsilon}_j = L^\top \varepsilon_j \sim N_n(0, I_n)$.
• Das Ergebnis ist eine de-korrelierte latente Datenmatrix $\tilde{Z}$, auf die beliebige Standard-DAG-Lernalgorithmen angewendet werden können.

C. Schätzverfahren (EM-Algorithmus und Pairwise MLE)

Da die latenten Variablen $Z$ nicht direkt beobachtbar sind (insbesondere bei diskreten Daten), wird ein zweiphasiger Ansatz verwendet:

1. Phase 1: Vor-Schätzung (Pre-estimation):

   • Schätzung der Schwellenwerte $T$ und Regressionskoeffizienten $\beta$ unter der Annahme $\Sigma = I_n$ (Unabhängigkeit), unter Verwendung eines Koordinatenanstiegs-Algorithmus.
   • Schätzung der Kovarianzmatrix $\Sigma$ mittels Pairwise Maximum Likelihood Estimation (MLE). Dabei wird die Korrelation zwischen Paaren von Einheiten basierend auf der gemeinsamen Verteilung aller Variablen (kontinuierlich und diskret) geschätzt. Es wird eine Block-diagonale Struktur für $\Sigma$ angenommen, um die Komplexität zu handhaben.

2. Phase 2: Latente Daten-Rekonstruktion und De-Korrelation:

   • Ein EM-Algorithmus wird verwendet, um die latenten Variablen $Z$ zu imputieren.
   • E-Schritt: Berechnung des Erwartungswerts der latenten Variablen gegeben die beobachteten diskreten Daten und die aktuellen Parameter. Da die direkte Integration hochdimensional ist, wird ein Gibbs-Sampler innerhalb der Blöcke der Kovarianzmatrix verwendet.
   • M-Schritt: Aktualisierung der Regressionskoeffizienten $\beta$ durch Regression der de-korrelierten latenten Daten auf die de-korrelierten Elternvariablen.
   • Abschließend wird die Cholesky-Transformation $L^\top$ auf die rekonstruierten latenten Daten angewendet, um $\tilde{Z}$ zu erhalten.

D. Struktur-Lernen

Auf den de-korrelierten Daten $\tilde{Z}$ werden Standard-Algorithmen (z. B. MMHC, PC, Copula-PC) angewendet. Um die Variabilität durch die Monte-Carlo-Schätzung im EM-Algorithmus zu reduzieren, werden zwei Strategien vorgeschlagen:

• Consensus: Mehrere de-korrelierte Datensätze werden generiert, und die resultierenden DAGs werden durch Mehrheitsvoting aggregiert.
• Average: Die de-korrelierten Datensätze werden gemittelt, und ein einziger DAG wird geschätzt.

3. Wichtige Beiträge

1. Erweiterung auf gemischte Daten: Das Paper generalisiert bestehende De-Korrelations-Frameworks (die nur für kontinuierliche oder nur für binäre Daten galten) auf eine Mischung aus kontinuierlichen und diskreten (multinomialen) Variablen.
2. Robustes Schätzverfahren: Entwicklung eines EM-Algorithmus mit Gibbs-Sampling zur Imputation latenter Variablen in Abhängigkeit von einer geschätzten Kovarianzstruktur zwischen Einheiten.
3. Anwendung auf scRNA-seq: Die Methode wird erfolgreich auf reale Daten zur Genregulationsnetzwerk-Inferenz (GRN) angewendet, wo Zellabhängigkeiten und gemischte Datentypen (z. B. diskretisierte Expressionslevel) eine große Rolle spielen.
4. Stabilitätsanalyse: Einführung eines Bootstrap-Verfahrens zur Quantifizierung der Unsicherheit der inferierten Kanten und zur Identifizierung hochkonfidenzierter regulatorischer Interaktionen.

4. Ergebnisse

Simulationen

• Kovarianzschätzung: Die Pairwise-MLE-Methode schätzt die Kovarianzmatrix $\Sigma$ robust, selbst wenn die Elternmengen der DAGs nur ungenau geschätzt wurden.
• DAG-Wiederherstellung: In Simulationen mit synthetischen und realen Benchmark-DAGs (aus dem bnlearn-Repository) übertraf das De-Korrelations-Framework (sowohl Consensus- als auch Average-Ansatz) die Baseline (direkte Anwendung auf abhängige Daten) signifikant.
  • Die Verbesserung war besonders stark bei hohen Variablenzahlen ($p=1000$) und bei Verwendung von Copula-PC, was zeigt, dass die Methode die Probleme der Diskretisierung und Abhängigkeit effektiv adressiert.
  • Der F1-Score für die Wiederherstellung der DAG-Struktur stieg in fast allen Szenarien an.

Anwendung auf Genregulationsnetzwerke (scRNA-seq)

• Daten: Analyse von 859 embryonalen Stammzellen und 51 wichtigen Genen (Chu et al., 2016).
• Modellgüte: Ein 10-faches Cross-Validation zeigte, dass das de-korrelierte Modell eine signifikant höhere Log-Likelihood auf Testdaten aufwies als das Baseline-Modell oder ein Modell ohne De-Korrelation. Dies bestätigt, dass die explizite Modellierung der Zellabhängigkeiten die Anpassungsgüte verbessert.
• Biologische Validität: Das inferierte Netzwerk enthielt hochkonfidente Kanten, die mit bekannten biologischen Interaktionen übereinstimmen (z. B. Interaktionen zwischen POU5F1, NANOG und ZFP42 sowie Feedback-Schleifen zwischen NODAL und LEFTY1). Viele dieser Kanten wurden durch Bootstrap-Analysen als stabil bestätigt.

5. Bedeutung und Fazit

Das Paper liefert einen wichtigen methodischen Fortschritt für die kausale Inferenz in modernen biologischen und sozialen Datensätzen.

• Paradigmenwechsel: Es zeigt, dass kausale Entdeckungsalgorithmen nicht neu entwickelt werden müssen, um mit abhängigen Daten umzugehen. Stattdessen kann das Problem durch eine vorgelagerte De-Korrelationstransformation gelöst werden.
• Biologische Relevanz: Für die Genomik ist dies entscheidend, da die Annahme der Unabhängigkeit von Zellen in scRNA-seq-Daten oft falsch ist. Die Methode ermöglicht präzisere Vorhersagen von Genregulationsnetzwerken, was für das Verständnis von Differenzierungsprozessen und Krankheitsmechanismen essenziell ist.
• Zukunftsausblick: Die Autoren sehen Potenzial in der Erweiterung auf nichtlineare SEMs und die Integration von Interventionsdaten, um die Identifizierbarkeit kausaler Strukturen weiter zu verbessern.

Zusammenfassend bietet das vorgeschlagene Framework eine robuste, skalierbare und biologisch plausible Lösung für das Problem der kausalen Entdeckung in komplexen, abhängigen und gemischten Datensätzen.