Clarify all-future model interpretation

This commit is contained in:
2026-06-27 10:53:55 +08:00
parent 6df0435f65
commit df6b45b95a
12 changed files with 393 additions and 4580 deletions

Binary file not shown.

View File

@@ -0,0 +1,393 @@
\documentclass[11pt]{ctexart}
\usepackage[margin=1in]{geometry}
\usepackage{amsmath, amssymb}
\usepackage{booktabs}
\usepackage{enumitem}
\usepackage{hyperref}
\title{Interpreting the All-Future DeepHealth Model\\
all-future DeepHealth 模型的解释边界}
\author{}
\date{}
\begin{document}
\maketitle
\begin{abstract}
The all-future DeepHealth model should be interpreted as a
history-conditioned incident disease and mortality risk model. At a query time
\(t\), the hidden state \(h(t)\) summarizes the observed history up to that
time, and the model estimates future first-occurrence risks for diseases in the
model vocabulary, together with future mortality risk. The model does not
directly estimate current clinical disease burden, organ damage, frailty,
disease severity, recurrence risk, or disease-specific weights. Those
quantities require labels, mappings, weights, or supervision that are not part
of the present model.
\medskip
all-future DeepHealth 模型应被解释为一个基于既往轨迹的未来新发疾病和死亡风险模型。
在查询时刻 \(t\),隐含状态 \(h(t)\) 汇总了截至该时刻的已观测历史;模型输出的是
模型词表内各疾病的未来首次发生风险,以及未来死亡风险。当前模型并不直接估计当前
临床疾病负担、器官损伤、衰弱程度、疾病严重度、复发风险或疾病特异性权重。这些量
都需要当前模型之外的标签、映射、权重或额外监督。
\end{abstract}
\section{English}
\subsection{Model object}
For individual \(i\), let \(\mathcal{H}_i(t)\) denote the observed history up to
query time \(t\). The all-future model produces
\[
h_i(t) = f_\theta(\mathcal{H}_i(t)),
\]
and, for each disease token \(d\) in the modeled disease vocabulary
\(\mathcal{D}\), estimates a future first-occurrence risk over a horizon
\(\tau\):
\[
p_{i,d}(t,\tau)
=
P_\theta\!\left(T_{i,d}\in(t,t+\tau]\mid h_i(t)\right),
\]
where \(T_{i,d}\) is the first observed occurrence time of disease \(d\). If the
model includes a death endpoint, it also estimates
\[
p_{i,\mathrm{death}}(t,\tau)
=
P_\theta\!\left(T_{i,\mathrm{death}}\in(t,t+\tau]\mid h_i(t)\right).
\]
The disease sequence is a first-occurrence sequence. Therefore, from the input
history itself we know
\[
m_{i,d}(t)=\mathbf{1}\{T_{i,d}\le t\}.
\]
This historical indicator is not learned by the model; it is read directly from
the event history.
\subsection{Masking already occurred diseases}
For a disease that has already occurred before or at \(t\), the model output for
that disease should not be interpreted as recurrence risk or current disease
activity. For future incident disease summaries, already occurred diseases
should be masked:
\[
p^{\mathrm{new}}_{i,d}(t,\tau)
=
\left[1-m_{i,d}(t)\right]p_{i,d}(t,\tau).
\]
\subsection{Directly supported model-derived quantities}
The current model directly supports the following quantities.
\paragraph{Disease-specific future first-occurrence risk.}
For each modeled disease \(d\),
\[
p^{\mathrm{new}}_{i,d}(t,\tau)
\]
is the estimated risk that disease \(d\) newly appears within the next
\(\tau\) years, conditional on the history summarized by \(h_i(t)\).
\paragraph{Future mortality risk.}
\[
p_{i,\mathrm{death}}(t,\tau)
\]
is the estimated probability of death within the next \(\tau\) years. Death is a
terminal endpoint and should not be treated as an ordinary disease burden
weight.
\paragraph{Probability of being alive with no new modeled disease.}
Using the model's disease-specific and death risks, one can summarize the
probability of no new modeled disease and survival over the next \(\tau\) years:
\[
S^{\mathrm{all}}_i(t,\tau)
=
\left[1-p_{i,\mathrm{death}}(t,\tau)\right]
\prod_{d\in\mathcal{D}}
\left[1-p^{\mathrm{new}}_{i,d}(t,\tau)\right].
\]
Equivalently, if the model is represented through cumulative hazards
\(\Lambda_{i,d}(t,\tau)=-\log[1-p_{i,d}(t,\tau)]\),
\[
S^{\mathrm{all}}_i(t,\tau)
=
\exp\left(
-\Lambda_{i,\mathrm{death}}(t,\tau)
-\sum_{d\in\mathcal{D}}
[1-m_{i,d}(t)]\Lambda_{i,d}(t,\tau)
\right).
\]
\paragraph{Probability of being alive with no new disease in a specified set.}
For any analyst-specified subset of disease tokens \(G\subseteq\mathcal{D}\),
\[
S^{G}_i(t,\tau)
=
\left[1-p_{i,\mathrm{death}}(t,\tau)\right]
\prod_{d\in G}
\left[1-p^{\mathrm{new}}_{i,d}(t,\tau)\right].
\]
This is a subset-level future disease-free survival summary. If \(G\) is called
an organ system, the disease-to-organ grouping is external to the model and
must not be described as a learned organ score.
\subsection{Historical counts are not model-derived burden}
One may count observed historical diseases:
\[
B^{\mathrm{history}}_i(t)
=
\sum_{d\in\mathcal{D}}m_{i,d}(t).
\]
This quantity is a direct count from the input sequence. It does not require the
model and should not be presented as a model-derived disease burden score.
Without disease severity labels or disease weights, it treats all disease
tokens equally.
\subsection{What the current model does not estimate}
The current model does not directly estimate:
\begin{itemize}[leftmargin=1.5em]
\item current clinical disease burden;
\item organ damage, organ age, or organ functional reserve;
\item frailty or frailty weights;
\item severity of a newly diagnosed disease;
\item recurrence risk after first occurrence;
\item relative clinical importance of different disease tokens.
\end{itemize}
These interpretations require additional labels, mappings, weights, or model
training objectives. Using the present all-future model to claim these
quantities would be over-interpretation.
\subsection{Post-onset prognosis for the same new disease}
For a disease \(d\) that newly occurs at time \(T_{i,d}\), the model cannot
infer the clinical severity of that disease itself. However, after the disease
has entered the history, one may query the model again and compare subsequent
future risks:
\[
p^{\mathrm{new}}_{i,e}(T_{i,d},\tau),\quad e\ne d,
\qquad
p_{i,\mathrm{death}}(T_{i,d},\tau).
\]
This supports a prognosis-oriented statement: the same incident diagnosis may
be followed by different future disease and mortality risk profiles in people
with different prior trajectories. It should not be described as direct disease
severity.
\subsection{Future extension with reliable recurrence data}
The above interpretation is constrained by the first-occurrence nature of the
current disease sequence. UK Biobank does not provide a reliable longitudinal
record of disease recurrence, relapse, repeated admissions, treatment
escalation, or episode-level severity for the modeled disease tokens. Therefore,
the present model cannot be used to estimate recurrence risk or ongoing disease
activity after first onset.
If reliable recurrence or repeated-event data were available, one could define a
different modeling target. Let \(N_{i,d}(t)\) be the counting process for all
episodes of disease \(d\), not only its first occurrence. A recurrence-aware
model could estimate the future increment
\[
P_\theta\!\left(N_{i,d}(t+\tau)-N_{i,d}(t)>0 \mid h_i(t)\right),
\]
or the expected number of future episodes
\[
E_\theta\!\left[N_{i,d}(t+\tau)-N_{i,d}(t)\mid h_i(t)\right].
\]
For individuals with \(m_{i,d}(t)=1\), this would support a genuine
post-onset interpretation of future recurrence or disease activity. For
individuals with \(m_{i,d}(t)=0\), it would remain an incident disease risk.
With such data, the model could also separate three clinically different
quantities:
\begin{itemize}[leftmargin=1.5em]
\item incident risk: the first future onset of a disease not yet observed;
\item recurrence or repeated-event risk: future episodes after a disease has
already occurred;
\item mortality risk: a terminal endpoint that competes with future
non-fatal events.
\end{itemize}
If recurrence episodes were linked to reliable episode-level severity labels,
such as hospitalization intensity, treatment escalation, or validated severity
grades, then a further supervised model could learn severity-aware prognosis.
These extensions would require new data and new training targets; they are not
available from the current all-future first-occurrence model.
\section{中文}
\subsection{模型对象}
对于个体 \(i\),令 \(\mathcal{H}_i(t)\) 表示查询时刻 \(t\) 之前已经观测到的历史。
all-future 模型产生
\[
h_i(t)=f_\theta(\mathcal{H}_i(t)),
\]
并且对模型疾病词表 \(\mathcal{D}\) 中的每个疾病 \(d\),估计未来 \(\tau\) 年内的首次
发生风险:
\[
p_{i,d}(t,\tau)
=
P_\theta\!\left(T_{i,d}\in(t,t+\tau]\mid h_i(t)\right),
\]
其中 \(T_{i,d}\) 是疾病 \(d\) 的首次观测发生时间。如果模型包含死亡终点,则同时估计
\[
p_{i,\mathrm{death}}(t,\tau)
=
P_\theta\!\left(T_{i,\mathrm{death}}\in(t,t+\tau]\mid h_i(t)\right).
\]
当前疾病序列是 first-occurrence 序列。因此,从输入历史本身即可得到
\[
m_{i,d}(t)=\mathbf{1}\{T_{i,d}\le t\}.
\]
这个历史发生指示量不是模型学出来的,而是直接从事件历史中读出的。
\subsection{已经发生疾病的 mask}
如果某个疾病在 \(t\) 之前或 \(t\) 时已经发生,那么该疾病对应的模型输出不应解释为复发
风险,也不应解释为当前疾病活跃程度。在汇总未来新发疾病风险时,应对已经发生过的疾病
进行 mask
\[
p^{\mathrm{new}}_{i,d}(t,\tau)
=
\left[1-m_{i,d}(t)\right]p_{i,d}(t,\tau).
\]
\subsection{当前模型直接支持的派生量}
当前模型直接支持以下几类量。
\paragraph{疾病层面的未来首次发生风险。}
对每个模型内疾病 \(d\)
\[
p^{\mathrm{new}}_{i,d}(t,\tau)
\]
表示在 \(h_i(t)\) 所总结的历史条件下,疾病 \(d\) 在未来 \(\tau\) 年内新发生的风险。
\paragraph{未来死亡风险。}
\[
p_{i,\mathrm{death}}(t,\tau)
\]
表示未来 \(\tau\) 年内死亡的概率。死亡是终末结局,不应作为普通疾病负担权重加入疾病
负担求和。
\paragraph{未来无新病且存活的概率。}
利用疾病层面风险和死亡风险,可以汇总未来 \(\tau\) 年内无任何模型内新病且存活的概率:
\[
S^{\mathrm{all}}_i(t,\tau)
=
\left[1-p_{i,\mathrm{death}}(t,\tau)\right]
\prod_{d\in\mathcal{D}}
\left[1-p^{\mathrm{new}}_{i,d}(t,\tau)\right].
\]
如果用累计 hazard 表示,令
\(\Lambda_{i,d}(t,\tau)=-\log[1-p_{i,d}(t,\tau)]\),则
\[
S^{\mathrm{all}}_i(t,\tau)
=
\exp\left(
-\Lambda_{i,\mathrm{death}}(t,\tau)
-\sum_{d\in\mathcal{D}}
[1-m_{i,d}(t)]\Lambda_{i,d}(t,\tau)
\right).
\]
\paragraph{未来无指定疾病集合新发且存活的概率。}
对于任意由分析者预先指定的疾病 token 集合 \(G\subseteq\mathcal{D}\),可以定义
\[
S^{G}_i(t,\tau)
=
\left[1-p_{i,\mathrm{death}}(t,\tau)\right]
\prod_{d\in G}
\left[1-p^{\mathrm{new}}_{i,d}(t,\tau)\right].
\]
这只是指定疾病集合层面的未来 disease-free survival 汇总。如果将 \(G\) 称为某个器官系统,
那么疾病到器官的分组来自模型外部,不能描述为模型学到的器官评分。
\subsection{历史计数不是模型派生的疾病负担}
可以计算历史已经发生过多少个模型内疾病:
\[
B^{\mathrm{history}}_i(t)
=
\sum_{d\in\mathcal{D}}m_{i,d}(t).
\]
但这个量只是从输入序列直接计数,不需要模型,因此不应表述为模型派生的疾病负担评分。
在没有疾病严重度标签或疾病权重时,它默认所有疾病 token 等价。
\subsection{当前模型不能估计什么}
当前模型不能直接估计:
\begin{itemize}[leftmargin=1.5em]
\item 当前临床疾病负担;
\item 器官损伤、器官年龄或器官功能储备;
\item 衰弱程度或衰弱权重;
\item 某个新发疾病本身的严重程度;
\item 首次发生后的复发风险;
\item 不同疾病 token 之间的相对临床重要性。
\end{itemize}
这些解释都需要额外标签、映射、权重或新的训练目标。用当前 all-future 模型直接声称这些量,
属于过分解读。
\subsection{同一新发疾病后的预后差异}
对于在 \(T_{i,d}\) 时刻新发生的疾病 \(d\),模型不能推断这个疾病本身的临床严重程度。
但是,当该疾病已经进入历史之后,可以再次查询模型,并比较之后的未来风险:
\[
p^{\mathrm{new}}_{i,e}(T_{i,d},\tau),\quad e\ne d,
\qquad
p_{i,\mathrm{death}}(T_{i,d},\tau).
\]
这支持一种预后层面的表述:同一个新发诊断出现在不同既往轨迹的人身上,可能对应不同的
后续新病和死亡风险结构。但这不应被表述为模型直接判断了疾病严重度。
\subsection{如果有可靠复发数据,未来可以做什么}
上述解释受到当前 first-occurrence 疾病序列的限制。UK Biobank 并没有为模型词表中的
疾病提供可靠的纵向复发、复燃、重复住院、治疗升级或 episode 层面严重程度记录。因此,
当前模型不能用于估计首次发病后的复发风险或持续疾病活动。
如果未来有可靠的复发或重复事件数据,可以定义一个不同的建模目标。令 \(N_{i,d}(t)\)
表示疾病 \(d\) 的所有事件计数过程,而不只是首次发生。一个考虑复发的模型可以估计未来
事件增量:
\[
P_\theta\!\left(N_{i,d}(t+\tau)-N_{i,d}(t)>0 \mid h_i(t)\right),
\]
或者估计未来事件次数的期望:
\[
E_\theta\!\left[N_{i,d}(t+\tau)-N_{i,d}(t)\mid h_i(t)\right].
\]
对于 \(m_{i,d}(t)=1\) 的个体,这才可以支持真正的发病后复发风险或疾病活动解释;对于
\(m_{i,d}(t)=0\) 的个体,它仍然对应未来首次发病风险。
在这样的数据条件下,模型可以区分三个临床上不同的量:
\begin{itemize}[leftmargin=1.5em]
\item 新发风险:尚未发生疾病的未来首次发生风险;
\item 复发或重复事件风险:疾病已经发生后的未来 episode 风险;
\item 死亡风险:与非致死事件竞争的终末结局风险。
\end{itemize}
如果复发 episode 还带有可靠的 episode 层面严重程度标签,例如住院强度、治疗升级或经过
验证的严重程度分级,那么还可以进一步训练有监督的严重程度相关预后模型。但这些扩展都
需要新的数据和新的训练目标,并不是当前 all-future first-occurrence 模型已经具备的能力。
\section{Recommended wording}
\paragraph{English.}
The all-future model is a history-conditioned incident disease and mortality
risk model. Its outputs support future disease-free survival summaries over the
modeled disease vocabulary, but do not directly quantify current disease
burden, organ damage, frailty, or disease severity.
\paragraph{中文。}
all-future 模型是基于既往轨迹的未来新发疾病和死亡风险模型。它的输出可以支持模型词表
范围内的未来无新病且存活概率汇总,但不能直接量化当前疾病负担、器官损伤、衰弱程度或
疾病严重度。
\end{document}

View File

@@ -1,157 +0,0 @@
from __future__ import annotations
import csv
import re
from pathlib import Path
LABELS = Path("labels.csv")
OUT = Path("organ_involvement_label_mapping.csv")
ORGANS = [
"brain_neurologic",
"heart",
"artery_vascular",
"immune",
"intestine_digestive",
"kidney",
"liver",
"lung",
"muscle_musculoskeletal",
"pancreas_endocrine",
"adipose_metabolic",
"female_reproductive",
"male_reproductive",
"neoplasm",
]
def _code_key(code: str) -> tuple[str, int, str]:
code = code.strip().upper()
match = re.match(r"^([A-Z])(\d{2})(?:\.?([A-Z0-9]+))?", code)
if not match:
raise ValueError(f"Invalid ICD-10 code: {code!r}")
letter, number, suffix = match.groups()
return letter, int(number), suffix or ""
def _in_range(code: str, start: str, end: str) -> bool:
c_letter, c_num, _ = _code_key(code)
s_letter, s_num, _ = _code_key(start)
e_letter, e_num, _ = _code_key(end)
if s_letter == e_letter:
return c_letter == s_letter and s_num <= c_num <= e_num
return (
(s_letter < c_letter < e_letter)
or (c_letter == s_letter and c_num >= s_num)
or (c_letter == e_letter and c_num <= e_num)
)
def _matches_any(code: str, ranges: list[tuple[str, str]]) -> bool:
return any(_in_range(code, start, end) for start, end in ranges)
def organ_for_icd10(code: str) -> tuple[str, str]:
code = code.strip().upper()
if not re.match(r"^[A-Z]\d{2}", code):
return "", "unmapped_non_icd10"
if _matches_any(code, [("C00", "D48")]):
return "neoplasm", "neoplasm_c00_d48"
if _matches_any(code, [("F00", "F09"), ("G00", "G99"), ("I60", "I69")]):
return "brain_neurologic", "f00_f09_g00_g99_i60_i69"
if _matches_any(code, [("I00", "I09"), ("I20", "I52")]):
return "heart", "i00_i09_i20_i52"
if _matches_any(code, [("I10", "I15"), ("I70", "I89"), ("I95", "I99")]):
return "artery_vascular", "i10_i15_i70_i89_i95_i99"
if _matches_any(code, [("A00", "B99"), ("D50", "D89")]):
return "immune", "a00_b99_d50_d89"
if _matches_any(code, [("J00", "J99")]):
return "lung", "j00_j99"
if _matches_any(code, [("K70", "K77")]):
return "liver", "k70_k77"
if _matches_any(code, [("K85", "K86"), ("E10", "E16")]):
return "pancreas_endocrine", "k85_k86_e10_e16"
if _matches_any(code, [("K00", "K69"), ("K78", "K84"), ("K87", "K93")]):
return "intestine_digestive", "k00_k69_k78_k84_k87_k93"
if _matches_any(code, [("N00", "N39")]):
return "kidney", "n00_n39"
if _matches_any(code, [("N70", "N98"), ("O00", "O99")]):
return "female_reproductive", "n70_n98_o00_o99"
if _matches_any(code, [("N40", "N53")]):
return "male_reproductive", "n40_n53"
if _matches_any(code, [("M00", "M99")]):
return "muscle_musculoskeletal", "m00_m99"
if _matches_any(code, [("E00", "E09"), ("E17", "E90")]):
return "adipose_metabolic", "e00_e09_e17_e90"
return "", "unmapped_no_organ_rule"
def main() -> None:
rows = []
with LABELS.open(encoding="utf-8") as f:
for i, line in enumerate(f):
line = line.strip()
if not line:
continue
parts = line.split(maxsplit=1)
code = parts[0].strip().upper()
name = parts[1].strip() if len(parts) > 1 else ""
organ_id, match_source = organ_for_icd10(code)
rows.append(
{
"token_id": i + 3,
"label_code": code,
"label_name": name,
"organ_id": organ_id,
"organ_label": organ_id,
"organ_weight": 1.0 if organ_id else 0.0,
"match_source": match_source,
"mapping_source": (
"organ-age-inspired clinical systems based on "
"Oh et al. Nature 2023; single-label ICD-10 rules"
),
}
)
fieldnames = [
"token_id",
"label_code",
"label_name",
"organ_id",
"organ_label",
"organ_weight",
"match_source",
"mapping_source",
]
with OUT.open("w", newline="", encoding="utf-8-sig") as f:
writer = csv.DictWriter(f, fieldnames=fieldnames)
writer.writeheader()
writer.writerows(rows)
mapped = [row for row in rows if row["organ_id"]]
print(f"labels: {len(rows)}")
print(f"mapped_labels: {len(mapped)}")
print(f"unmapped_labels: {len(rows) - len(mapped)}")
print(f"organs: {', '.join(ORGANS)}")
print(f"output: {OUT}")
if __name__ == "__main__":
main()

View File

@@ -1,215 +0,0 @@
from __future__ import annotations
import csv
from pathlib import Path
LABELS = Path("labels.csv")
OUT = Path("uk_hfrs_label_mapping.csv")
MISSING_OUT = Path("uk_hfrs_missing_label_codes.csv")
# Source: Gilbert T, Neuburger J, Kraindler J, et al. Development and
# validation of a Hospital Frailty Risk Score focusing on older people in
# acute care settings using electronic hospital records. Lancet. 2018.
# Supplementary appendix, Table A2.
UK_HFRS_WEIGHTS = {
"F00": 7.1,
"G81": 4.4,
"G30": 4.0,
"I69": 3.7,
"R29": 3.6,
"N39": 3.2,
"F05": 3.2,
"W19": 3.2,
"S00": 3.2,
"R31": 3.0,
"B96": 2.9,
"R41": 2.7,
"R26": 2.6,
"I67": 2.6,
"R56": 2.6,
"R40": 2.5,
"T83": 2.4,
"S06": 2.4,
"S42": 2.3,
"E87": 2.3,
"M25": 2.3,
"E86": 2.3,
"R54": 2.2,
"Z50": 2.1,
"F03": 2.1,
"W18": 2.1,
"Z75": 2.0,
"F01": 2.0,
"S80": 2.0,
"L03": 2.0,
"H54": 1.9,
"E53": 1.9,
"Z60": 1.8,
"G20": 1.8,
"R55": 1.8,
"S22": 1.8,
"K59": 1.8,
"N17": 1.8,
"L89": 1.7,
"Z22": 1.7,
"B95": 1.7,
"L97": 1.6,
"R44": 1.6,
"K26": 1.6,
"I95": 1.6,
"N19": 1.6,
"A41": 1.6,
"Z87": 1.5,
"J96": 1.5,
"X59": 1.5,
"M19": 1.5,
"G40": 1.5,
"M81": 1.4,
"S72": 1.4,
"S32": 1.4,
"E16": 1.4,
"R94": 1.4,
"N18": 1.4,
"R33": 1.3,
"R69": 1.3,
"N28": 1.3,
"R32": 1.2,
"G31": 1.2,
"Y95": 1.2,
"S09": 1.2,
"R45": 1.2,
"G45": 1.2,
"Z74": 1.1,
"M79": 1.1,
"W06": 1.1,
"S01": 1.1,
"A04": 1.1,
"A09": 1.1,
"J18": 1.1,
"J69": 1.0,
"R47": 1.0,
"E55": 1.0,
"Z93": 1.0,
"R02": 1.0,
"R63": 0.9,
"H91": 0.9,
"W10": 0.9,
"W01": 0.9,
"E05": 0.9,
"M41": 0.9,
"R13": 0.8,
"Z99": 0.8,
"U80": 0.8,
"M80": 0.8,
"K92": 0.8,
"I63": 0.8,
"N20": 0.7,
"F10": 0.7,
"Y84": 0.7,
"R00": 0.7,
"J22": 0.7,
"Z73": 0.6,
"R79": 0.6,
"Z91": 0.5,
"S51": 0.5,
"F32": 0.5,
"M48": 0.5,
"E83": 0.4,
"M15": 0.4,
"D64": 0.4,
"L08": 0.4,
"R11": 0.3,
"K52": 0.3,
"R50": 0.1,
}
def _read_labels(path: Path) -> list[dict[str, str | int]]:
rows: list[dict[str, str | int]] = []
with path.open(encoding="utf-8") as f:
for i, line in enumerate(f):
line = line.strip()
if not line:
continue
parts = line.split(maxsplit=1)
code = parts[0].strip().upper()
name = parts[1].strip() if len(parts) > 1 else ""
rows.append({"token_id": i + 3, "label_code": code, "label_name": name})
return rows
def main() -> None:
labels = _read_labels(LABELS)
label_codes = {str(row["label_code"]) for row in labels}
missing = sorted(set(UK_HFRS_WEIGHTS) - label_codes)
rows = []
for row in labels:
code = str(row["label_code"])
weight = float(UK_HFRS_WEIGHTS.get(code, 0.0))
rows.append(
{
**row,
"hfrs_dimension_id": "hfrs_weighted_disease_expression",
"hfrs_dimension": "DeepHealth HFRS-weighted disease expression",
"hfrs_key_area": "UK-HFRS",
"hfrs_weight": weight,
"hfrs_source": (
"Gilbert et al. Lancet 2018 supplementary appendix Table A2"
),
"match_source": "exact_three_character_icd10" if weight else "not_in_hfrs",
}
)
fieldnames = [
"token_id",
"label_code",
"label_name",
"hfrs_dimension_id",
"hfrs_dimension",
"hfrs_key_area",
"hfrs_weight",
"hfrs_source",
"match_source",
]
with OUT.open("w", newline="", encoding="utf-8-sig") as f:
writer = csv.DictWriter(f, fieldnames=fieldnames)
writer.writeheader()
writer.writerows(rows)
with MISSING_OUT.open("w", newline="", encoding="utf-8-sig") as f:
writer = csv.DictWriter(
f,
fieldnames=[
"hfrs_source_code",
"hfrs_weight",
"missing_reason",
"hfrs_source",
],
)
writer.writeheader()
for code in missing:
writer.writerow(
{
"hfrs_source_code": code,
"hfrs_weight": UK_HFRS_WEIGHTS[code],
"missing_reason": "not_present_in_labels_csv",
"hfrs_source": (
"Gilbert et al. Lancet 2018 supplementary appendix Table A2"
),
}
)
nonzero = sum(1 for row in rows if float(row["hfrs_weight"]) != 0.0)
print(f"labels: {len(rows)}")
print(f"uk_hfrs_codes: {len(UK_HFRS_WEIGHTS)}")
print(f"matched_nonzero_labels: {nonzero}")
print(f"missing_hfrs_codes: {len(missing)}")
print(f"output: {OUT}")
print(f"missing_output: {MISSING_OUT}")
if __name__ == "__main__":
main()

View File

@@ -1,407 +0,0 @@
from __future__ import annotations
from dataclasses import dataclass
from pathlib import Path
from typing import Any, Sequence
import numpy as np
import torch
import torch.nn.functional as F
from eval_data import load_sequence_eval_dataset
from evaluate_auc_v2 import (
build_model_from_dataset,
load_checkpoint_state_dict,
load_json_config,
load_model_state,
resolve_dist_mode_for_checkpoint,
validate_dataset_metadata,
)
from targets import CHECKUP_IDX, NO_EVENT_IDX, PAD_IDX
@dataclass(frozen=True)
class DeepHealthContext:
model: torch.nn.Module
dataset: Any
cfg: dict[str, Any]
dist_mode: str
device: torch.device
run_path: Path
@dataclass(frozen=True)
class DiseaseExpressionResult:
disease_ids: np.ndarray
expression: np.ndarray
t_query: float
@dataclass(frozen=True)
class OrganInvolvementResult:
organ_ids: list[str]
involvement: np.ndarray
disease_ids: np.ndarray
expression: np.ndarray
t_query: float
@dataclass(frozen=True)
class FrailtyRiskResult:
frailty_risk_index: float
disease_ids: np.ndarray
expression: np.ndarray
weights: np.ndarray
t_query: float
def load_deephealth_context(
run_path: str | Path,
*,
device: str | torch.device | None = None,
) -> DeepHealthContext:
run_path = Path(run_path)
config_path = run_path / "train_config.json"
model_ckpt_path = run_path / "best_model.pt"
if not config_path.exists():
raise FileNotFoundError(f"train_config.json not found in {run_path}")
if not model_ckpt_path.exists():
raise FileNotFoundError(f"best_model.pt not found in {run_path}")
cfg = load_json_config(config_path)
model_target_mode = str(cfg.get("model_target_mode", "next_token")).lower()
if model_target_mode == "next_token":
raise RuntimeError(
"Disease expression indices require an all_future checkpoint because "
"they use p_d(h, Delta). The provided run is model_target_mode='next_token'."
)
if model_target_mode != "all_future":
raise ValueError(
"train_config.json model_target_mode must be all_future, got "
f"{model_target_mode!r}."
)
device_obj = torch.device(
device if device is not None else ("cuda" if torch.cuda.is_available() else "cpu")
)
if device_obj.type == "cuda" and not torch.cuda.is_available():
raise RuntimeError(f"Requested device {device_obj}, but CUDA is not available.")
dataset = load_sequence_eval_dataset(
model_target_mode="all_future",
data_prefix=cfg.get("data_prefix", "ukb"),
labels_file=cfg.get("labels_file", "labels.csv"),
no_event_interval_years=float(cfg.get("no_event_interval_years", 5.0)),
include_no_event_in_uts_target=bool(
cfg.get("include_no_event_in_uts_target", False)
),
min_history_events=int(cfg.get("all_future_min_history_events", 1)),
min_future_events=int(cfg.get("all_future_min_future_events", 1)),
extra_info_types=cfg.get("extra_info_types", None),
)
validate_dataset_metadata(dataset, cfg)
state_dict = load_checkpoint_state_dict(model_ckpt_path, map_location="cpu")
dist_mode = resolve_dist_mode_for_checkpoint(
str(cfg.get("dist_mode", "exponential")),
state_dict,
)
if dist_mode not in {"exponential", "weibull", "mixed"}:
raise ValueError(f"Unsupported dist_mode={dist_mode!r}")
cfg_model = dict(cfg)
cfg_model["dist_mode"] = dist_mode
model = build_model_from_dataset(_ConfigNamespace(), cfg_model, dataset)
load_model_state(model, state_dict)
model.eval().to(device_obj)
return DeepHealthContext(
model=model,
dataset=dataset,
cfg=cfg,
dist_mode=dist_mode,
device=device_obj,
run_path=run_path,
)
@torch.inference_mode()
def compute_disease_expression(
*,
run_path: str | Path,
disease_ids: Sequence[int] | np.ndarray,
event_seq: Sequence[int] | np.ndarray,
time_seq: Sequence[float] | np.ndarray,
sex: int,
other_type: Sequence[int] | np.ndarray,
other_value: Sequence[float] | np.ndarray,
other_value_kind: Sequence[int] | np.ndarray,
other_time: Sequence[float] | np.ndarray,
t_query: float,
device: str | torch.device | None = None,
context: DeepHealthContext | None = None,
) -> DiseaseExpressionResult:
ctx = context or load_deephealth_context(run_path, device=device)
disease_ids_arr = np.asarray(disease_ids, dtype=np.int64)
expression = model_implied_disease_expression(
ctx=ctx,
disease_ids=disease_ids_arr,
event_seq=event_seq,
time_seq=time_seq,
sex=sex,
other_type=other_type,
other_value=other_value,
other_value_kind=other_value_kind,
other_time=other_time,
t_query=float(t_query),
)
return DiseaseExpressionResult(
disease_ids=disease_ids_arr.copy(),
expression=expression,
t_query=float(t_query),
)
def compute_organ_involvement_from_expression(
*,
expression: Sequence[float] | np.ndarray,
organ_matrix: np.ndarray,
) -> np.ndarray:
z = np.asarray(expression, dtype=np.float64)
A = np.asarray(organ_matrix, dtype=np.float64)
if A.ndim != 2:
raise ValueError(f"organ_matrix must be 2D, got shape {A.shape}")
if z.ndim != 1 or A.shape[1] != z.size:
raise ValueError(
"expression must be 1D and match organ_matrix columns, got "
f"{z.shape} and {A.shape}"
)
intensity = -np.log1p(-np.clip(z, 0.0, 1.0 - 1e-7))
return -np.expm1(-(A @ intensity))
def compute_frailty_risk_from_expression(
*,
expression: Sequence[float] | np.ndarray,
hfrs_weights: Sequence[float] | np.ndarray,
) -> float:
z = np.asarray(expression, dtype=np.float64)
w = np.asarray(hfrs_weights, dtype=np.float64)
if z.shape != w.shape:
raise ValueError(f"expression and hfrs_weights shape mismatch: {z.shape} vs {w.shape}")
return float(np.dot(w, z))
@torch.inference_mode()
def model_implied_disease_expression(
*,
ctx: DeepHealthContext,
disease_ids: np.ndarray,
event_seq: Sequence[int] | np.ndarray,
time_seq: Sequence[float] | np.ndarray,
sex: int,
other_type: Sequence[int] | np.ndarray,
other_value: Sequence[float] | np.ndarray,
other_value_kind: Sequence[int] | np.ndarray,
other_time: Sequence[float] | np.ndarray,
t_query: float,
) -> np.ndarray:
disease_ids = np.asarray(disease_ids, dtype=np.int64)
_validate_disease_ids(ctx, disease_ids)
event_seq_arr, time_seq_arr = _validate_event_inputs(event_seq, time_seq)
other_type_arr, other_value_arr, other_value_kind_arr, other_time_arr = (
_validate_other_inputs(other_type, other_value, other_value_kind, other_time)
)
grid = build_readout_grid(
event_seq=event_seq_arr,
time_seq=time_seq_arr,
other_type=other_type_arr,
other_time=other_time_arr,
t_query=float(t_query),
)
if grid.size == 0:
return np.zeros(disease_ids.size, dtype=np.float64)
end_times = np.concatenate([grid[1:], np.asarray([t_query], dtype=np.float32)])
deltas = np.maximum(end_times - grid, 0.0).astype(np.float32)
valid = deltas > 0
if not np.any(valid):
return np.zeros(disease_ids.size, dtype=np.float64)
hidden = query_hidden(
ctx=ctx,
event_seq=event_seq_arr,
time_seq=time_seq_arr,
sex=sex,
other_type=other_type_arr,
other_value=other_value_arr,
other_value_kind=other_value_kind_arr,
other_time=other_time_arr,
query_times=grid[valid].astype(np.float32),
)
interval_prob = probabilities_from_hidden(
ctx=ctx,
hidden=hidden,
disease_ids=disease_ids,
deltas=deltas[valid],
)
survival = np.prod(1.0 - np.clip(interval_prob, 0.0, 1.0), axis=0)
return (1.0 - survival).astype(np.float64, copy=False)
def build_readout_grid(
*,
event_seq: np.ndarray,
time_seq: np.ndarray,
other_type: np.ndarray,
other_time: np.ndarray,
t_query: float,
) -> np.ndarray:
event_mask = (event_seq > PAD_IDX) & (time_seq <= np.float32(t_query))
other_mask = (other_type > 0) & (other_time <= np.float32(t_query))
times = np.concatenate(
[
time_seq[event_mask].astype(np.float32, copy=False),
other_time[other_mask].astype(np.float32, copy=False),
]
)
if times.size == 0:
return np.zeros(0, dtype=np.float32)
return np.unique(times)
@torch.inference_mode()
def query_hidden(
*,
ctx: DeepHealthContext,
event_seq: np.ndarray,
time_seq: np.ndarray,
sex: int,
other_type: np.ndarray,
other_value: np.ndarray,
other_value_kind: np.ndarray,
other_time: np.ndarray,
query_times: np.ndarray,
) -> torch.Tensor:
if query_times.ndim != 1:
raise ValueError("query_times must be 1D.")
batch_size = int(query_times.size)
if batch_size == 0:
return torch.empty(0, ctx.model.n_embd, device=ctx.device)
event = torch.from_numpy(event_seq[None, :].repeat(batch_size, axis=0)).long()
times = torch.from_numpy(time_seq[None, :].repeat(batch_size, axis=0)).float()
other_t = torch.from_numpy(other_type[None, :].repeat(batch_size, axis=0)).long()
other_v = torch.from_numpy(other_value[None, :].repeat(batch_size, axis=0)).float()
other_k = torch.from_numpy(
other_value_kind[None, :].repeat(batch_size, axis=0)
).long()
other_tm = torch.from_numpy(other_time[None, :].repeat(batch_size, axis=0)).float()
sex_t = torch.full((batch_size,), int(sex), dtype=torch.long)
tq = torch.from_numpy(query_times.astype(np.float32, copy=False)).float()
event = event.to(ctx.device)
return ctx.model(
event_seq=event,
time_seq=times.to(ctx.device),
sex=sex_t.to(ctx.device),
padding_mask=event > PAD_IDX,
t_query=tq.to(ctx.device),
other_type=other_t.to(ctx.device),
other_value=other_v.to(ctx.device),
other_value_kind=other_k.to(ctx.device),
other_time=other_tm.to(ctx.device),
target_mode="all_future",
)
@torch.inference_mode()
def probabilities_from_hidden(
*,
ctx: DeepHealthContext,
hidden: torch.Tensor,
disease_ids: np.ndarray,
deltas: np.ndarray,
) -> np.ndarray:
if hidden.ndim != 2:
raise ValueError(f"hidden must have shape (N, H), got {tuple(hidden.shape)}")
if deltas.ndim != 1 or deltas.size != hidden.shape[0]:
raise ValueError(
"deltas must be 1D with the same length as hidden rows, got "
f"{deltas.shape} vs {tuple(hidden.shape)}"
)
ids = torch.as_tensor(disease_ids, dtype=torch.long, device=ctx.device)
logits = ctx.model.calc_risk(hidden)[:, ids]
rate = F.softplus(logits).clamp_min(1e-8)
delta_t = torch.as_tensor(deltas, dtype=rate.dtype, device=ctx.device).clamp_min(0)
if ctx.dist_mode == "weibull":
rho = ctx.model.calc_weibull_rho(hidden)[:, ids]
exposure = torch.pow(delta_t[:, None], rho)
elif ctx.dist_mode == "mixed":
exposure = delta_t[:, None].expand_as(rate)
death_idx = int(getattr(ctx.model, "death_idx", getattr(ctx.model, "vocab_size", 0) - 1))
death_cols = [j for j, token in enumerate(disease_ids.tolist()) if int(token) == death_idx]
if death_cols:
death_rho = ctx.model.calc_death_rho(hidden)
for col in death_cols:
exposure[:, int(col)] = torch.pow(delta_t, death_rho)
else:
exposure = delta_t[:, None].expand_as(rate)
prob = -torch.expm1(-rate * exposure)
return prob.detach().cpu().numpy().astype(np.float64, copy=False)
class _ConfigNamespace:
def __getattr__(self, _name: str) -> None:
return None
def _validate_disease_ids(ctx: DeepHealthContext, disease_ids: np.ndarray) -> None:
if disease_ids.ndim != 1 or disease_ids.size == 0:
raise ValueError("disease_ids must be a non-empty 1D array.")
vocab_size = int(getattr(ctx.model, "vocab_size", ctx.model.risk_head.out_features))
if np.any(disease_ids < 0) or np.any(disease_ids >= vocab_size):
raise ValueError(f"disease_ids must be in [0, {vocab_size}), got {disease_ids}")
def _validate_event_inputs(
event_seq: Sequence[int] | np.ndarray,
time_seq: Sequence[float] | np.ndarray,
) -> tuple[np.ndarray, np.ndarray]:
events = np.asarray(event_seq, dtype=np.int64)
times = np.asarray(time_seq, dtype=np.float32)
if events.ndim != 1 or times.ndim != 1:
raise ValueError("event_seq and time_seq must be 1D.")
if events.shape != times.shape:
raise ValueError(
f"event_seq and time_seq must have the same shape, got {events.shape} vs {times.shape}"
)
if events.size == 0:
raise ValueError("event_seq must contain at least one token.")
if not np.all(np.isfinite(times)):
raise ValueError("time_seq contains non-finite values.")
return events, times
def _validate_other_inputs(
other_type: Sequence[int] | np.ndarray,
other_value: Sequence[float] | np.ndarray,
other_value_kind: Sequence[int] | np.ndarray,
other_time: Sequence[float] | np.ndarray,
) -> tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]:
typ = np.asarray(other_type, dtype=np.int64)
val = np.asarray(other_value, dtype=np.float32)
kind = np.asarray(other_value_kind, dtype=np.int64)
tm = np.asarray(other_time, dtype=np.float32)
if not (typ.shape == val.shape == kind.shape == tm.shape):
raise ValueError(
"other_type, other_value, other_value_kind, and other_time must "
f"have the same shape, got {typ.shape}, {val.shape}, {kind.shape}, {tm.shape}."
)
if typ.ndim != 1:
raise ValueError("other_* inputs must be 1D.")
if not np.all(np.isfinite(tm)):
raise ValueError("other_time contains non-finite values.")
return typ, val, kind, tm

View File

@@ -1,161 +0,0 @@
\documentclass[11pt]{article}
\usepackage[margin=1in]{geometry}
\usepackage{amsmath, amssymb}
\usepackage{booktabs}
\usepackage{enumitem}
\usepackage{hyperref}
\title{DeepHealth Disease Expression, Organ Involvement, and Frailty Risk Indices}
\author{}
\date{}
\begin{document}
\maketitle
\begin{abstract}
DeepHealth provides a query-time hidden state \(h(t)\) and disease-specific
risk functions \(p_d(h,\Delta)\). We use these outputs to define a continuous
disease expression rate \(z_d(t)\). This quantity should be interpreted as how
much disease \(d\) is model-implied to have formed or expressed by query time
\(t\), not as true physiological damage. Based on \(z_d(t)\), we define two
downstream indices: an organ involvement index, which summarizes whether an
organ-age-inspired clinical system is involved by any related disease process,
and a DeepHealth-HFRS frailty risk index, which is the original UK-HFRS weighted
sum with binary disease occurrence replaced by continuous disease expression.
\end{abstract}
\section{Disease Expression Rate}
For a patient queried at time \(t\), let the historical readout times be
\[
t_0 < t_1 < \cdots < t_n \le t,\qquad t_{n+1}=t.
\]
For each interval \([t_i,t_{i+1}]\), DeepHealth produces a hidden state
\(h_i=h(t_i)\) and an interval risk
\[
q_{d,i}(t)=p_d(h_i,t_{i+1}-t_i).
\]
The model-implied disease expression rate is defined by noisy-or accumulation:
\[
z_d(t)
=
1-\prod_{i=0}^{n}\left[1-q_{d,i}(t)\right].
\]
Informally, \(z_d(t)\) is the degree to which disease \(d\) is expressed in the
patient by time \(t\). Unlike a raw diagnosis indicator, it is continuous and
can reflect heterogeneity within the same ICD label.
\section{Organ Involvement Index}
The organ index is not a frailty score, health reserve score, or organ age. It
is an organ involvement index. Let \(\mathcal{D}_k\) be the set of diseases
assigned to organ/system \(k\). Define disease expression intensity as
\[
\Lambda_d(t)=-\log\left[1-z_d(t)\right].
\]
The equal-weight organ involvement index is
\[
O_k(t)
=
1-\exp\left(
-\sum_{d\in\mathcal{D}_k}\Lambda_d(t)
\right).
\]
Equivalently,
\[
O_k(t)
=
1-
\prod_{d\in\mathcal{D}_k}
\left[1-z_d(t)\right].
\]
Thus \(O_k(t)\in[0,1]\) is the probability-like degree to which organ/system
\(k\) is involved by at least one related disease process. In the current
version all diseases assigned to the same organ are equally weighted; this is a
first-stage structural definition. Future versions can introduce
organ-specific disease weights \(\alpha_{k,d}\):
\[
O_k(t)
=
1-\exp\left(
-\sum_{d\in\mathcal{D}_k}\alpha_{k,d}\Lambda_d(t)
\right).
\]
\section{Organ List}
The organ/system categories are inspired by organ-age studies, especially
organ-specific plasma proteomic aging models, and are adapted to ICD disease
labels. The current list is:
\begin{center}
\begin{tabular}{ll}
\toprule
ID & Label \\
\midrule
brain\_neurologic & Brain and neurologic system \\
heart & Heart \\
artery\_vascular & Artery and vascular system \\
immune & Immune and infection-related system \\
intestine\_digestive & Intestine and digestive system \\
kidney & Kidney and urinary system \\
liver & Liver \\
lung & Lung and respiratory system \\
muscle\_musculoskeletal & Muscle and musculoskeletal system \\
pancreas\_endocrine & Pancreas and endocrine system \\
adipose\_metabolic & Adipose and metabolic system \\
female\_reproductive & Female reproductive system \\
male\_reproductive & Male reproductive system \\
neoplasm & Neoplasm \\
\bottomrule
\end{tabular}
\end{center}
The neoplasm category is retained as a disease-system category rather than
forced into a single anatomical organ. Sex-specific reproductive diseases are
separated into female and male reproductive systems.
\section{DeepHealth-HFRS Frailty Risk Index}
The original UK-HFRS is a weighted sum over binary disease occurrence:
\[
\operatorname{HFRS}^{\mathrm{obs}}(t)
=
\sum_{d\in\mathcal{D}_{\mathrm{HFRS}}}
w^{\mathrm{HFRS}}_d\,o_d(t),
\qquad
o_d(t)\in\{0,1\}.
\]
DeepHealth-HFRS keeps the published UK-HFRS weights and replaces the binary
disease state with the continuous DeepHealth disease expression rate:
\[
\operatorname{HFRS}^{\mathrm{DH}}(t)
=
\sum_{d\in\mathcal{D}_{\mathrm{HFRS}}}
w^{\mathrm{HFRS}}_d\,z_d(t),
\qquad
z_d(t)\in[0,1].
\]
This is a natural continuous extension of the original HFRS, so it can still be
called a frailty risk index. The semantic change is not the HFRS weight system;
the change is the disease state variable.
\section{Current Implementation}
The current code computes historical current-state indices only. No future
horizon is used. For each landmark age \(t\), it outputs:
\begin{itemize}[leftmargin=*]
\item \(z_d(t)\) internally as model-implied disease expression;
\item \(O_k(t)\) as equal-weight organ involvement;
\item \(\operatorname{HFRS}^{\mathrm{DH}}(t)\) as DeepHealth-HFRS frailty
risk.
\end{itemize}
The output table uses the columns
\[
\texttt{index\_type},\quad
\texttt{index\_id},\quad
\texttt{index\_label},\quad
\texttt{index\_value}.
\]
\end{document}

View File

@@ -1,142 +0,0 @@
\documentclass[11pt]{ctexart}
\usepackage[margin=1in]{geometry}
\usepackage{amsmath, amssymb}
\usepackage{booktabs}
\usepackage{enumitem}
\usepackage{hyperref}
\title{DeepHealth 疾病表达率、器官受累指数与衰弱风险指数}
\author{}
\date{}
\begin{document}
\maketitle
\begin{abstract}
DeepHealth 在查询时刻 \(t\) 输出隐含状态 \(h(t)\),并给出疾病风险函数
\(p_d(h,\Delta)\)。我们首先定义连续的疾病表达率 \(z_d(t)\):它表示模型认为疾病 \(d\)
截至 \(t\) 在该个体身上形成或表达了多少,而不是疾病造成的真实损害。基于 \(z_d(t)\)
本文定义两类指数:器官受累指数和 DeepHealth-HFRS 衰弱风险指数。前者表示器官/系统是否被相关疾病过程累及;
后者是原版 UK-HFRS 的自然连续化,即用连续疾病表达率替代二值疾病发生状态。
\end{abstract}
\section{疾病表达率}
设历史 readout 时间为
\[
t_0 < t_1 < \cdots < t_n \le t,\qquad t_{n+1}=t.
\]
在区间 \([t_i,t_{i+1}]\) 上,模型给出疾病 \(d\) 的区间风险
\[
q_{d,i}(t)=p_d(h(t_i),t_{i+1}-t_i).
\]
疾病表达率定义为
\[
z_d(t)
=
1-\prod_{i=0}^{n}\left[1-q_{d,i}(t)\right].
\]
直观上,\(z_d(t)\) 表示“这个病在该个体身上形成或表达了多少”。它不是二值诊断记录,
因此可以表达同一 ICD 标签下的个体异质性。
\section{器官受累指数}
器官指数不定义为器官年龄、器官健康储备或器官衰弱,而定义为器官受累指数。设 \(\mathcal{D}_k\)
是归属于器官/系统 \(k\) 的疾病集合。定义疾病表达强度
\[
\Lambda_d(t)=-\log\left[1-z_d(t)\right].
\]
当前版本使用等权器官受累定义:
\[
O_k(t)
=
1-\exp\left(
-\sum_{d\in\mathcal{D}_k}\Lambda_d(t)
\right),
\]
等价于
\[
O_k(t)
=
1-
\prod_{d\in\mathcal{D}_k}
\left[1-z_d(t)\right].
\]
因此 \(O_k(t)\in[0,1]\),表示器官/系统 \(k\) 被至少一个相关疾病过程累及的概率型程度。
当前所有疾病在同一器官内等权;后续可扩展为带疾病权重的形式:
\[
O_k(t)
=
1-\exp\left(
-\sum_{d\in\mathcal{D}_k}\alpha_{k,d}\Lambda_d(t)
\right).
\]
\section{器官列表}
当前器官/系统列表参考器官年龄研究中的 organ-age-inspired systems并结合 ICD 疾病标签空间调整:
\begin{center}
\begin{tabular}{ll}
\toprule
ID & 含义 \\
\midrule
brain\_neurologic & 脑与神经系统 \\
heart & 心脏 \\
artery\_vascular & 动脉与血管系统 \\
immune & 免疫与感染相关系统 \\
intestine\_digestive & 肠道与消化系统 \\
kidney & 肾脏与泌尿系统 \\
liver & 肝脏 \\
lung & 肺与呼吸系统 \\
muscle\_musculoskeletal & 肌肉骨骼系统 \\
pancreas\_endocrine & 胰腺与内分泌系统 \\
adipose\_metabolic & 脂肪与代谢系统 \\
female\_reproductive & 女性生殖系统 \\
male\_reproductive & 男性生殖系统 \\
neoplasm & 肿瘤 \\
\bottomrule
\end{tabular}
\end{center}
肿瘤作为疾病系统单独保留,不强行归入某个单一器官。男女生殖系统单独拆分。
\section{DeepHealth-HFRS 衰弱风险指数}
原版 UK-HFRS 是二值疾病发生状态的加权和:
\[
\operatorname{HFRS}^{\mathrm{obs}}(t)
=
\sum_{d\in\mathcal{D}_{\mathrm{HFRS}}}
w^{\mathrm{HFRS}}_d o_d(t),
\qquad
o_d(t)\in\{0,1\}.
\]
DeepHealth-HFRS 保留原版 UK-HFRS 权重,只把疾病状态从二值观测替换为连续疾病表达率:
\[
\operatorname{HFRS}^{\mathrm{DH}}(t)
=
\sum_{d\in\mathcal{D}_{\mathrm{HFRS}}}
w^{\mathrm{HFRS}}_d z_d(t),
\qquad
z_d(t)\in[0,1].
\]
因此 DeepHealth-HFRS 仍然可以称为衰弱风险指数;它是原版 HFRS 的自然连续化。
\section{当前实现}
当前代码只计算历史当前状态,不再使用未来 horizon。每个 landmark age \(t\) 输出:
\begin{itemize}[leftmargin=*]
\item 内部疾病表达率 \(z_d(t)\)
\item 等权器官受累指数 \(O_k(t)\)
\item DeepHealth-HFRS 衰弱风险指数 \(\operatorname{HFRS}^{\mathrm{DH}}(t)\)
\end{itemize}
输出表使用
\[
\texttt{index\_type},\quad
\texttt{index\_id},\quad
\texttt{index\_label},\quad
\texttt{index\_value}.
\]
\end{document}

View File

@@ -1,717 +0,0 @@
from __future__ import annotations
import argparse
import multiprocessing as mp
import time
from concurrent.futures import ProcessPoolExecutor, as_completed
from pathlib import Path
from typing import Any, Iterable
import numpy as np
import pandas as pd
import torch
from torch.nn.utils.rnn import pad_sequence
from torch.utils.data import DataLoader, IterableDataset, get_worker_info
from tqdm.auto import tqdm
from burden_index import (
build_readout_grid,
load_deephealth_context,
probabilities_from_hidden,
)
from evaluate_auc_v2 import make_eval_indices, parse_float_list
from targets import CHECKUP_IDX, NO_EVENT_IDX, PAD_IDX
def _parse_landmark_ages(args: argparse.Namespace) -> np.ndarray:
explicit = parse_float_list(args.landmark_ages)
if explicit:
ages = np.asarray(explicit, dtype=np.float32)
else:
ages = np.arange(
float(args.landmark_start),
float(args.landmark_stop) + 1e-6,
float(args.landmark_step),
dtype=np.float32,
)
if ages.size == 0:
raise ValueError("No landmark ages were provided.")
return ages
def _parse_devices(args: argparse.Namespace) -> list[str | None]:
if args.devices is not None and str(args.devices).strip():
devices = [x.strip() for x in str(args.devices).split(",") if x.strip()]
if not devices:
raise ValueError("--devices was provided but no devices were parsed.")
return devices
return [args.device]
def _load_index_matrices(
*,
organ_mapping_csv: Path,
hfrs_mapping_csv: Path,
) -> tuple[np.ndarray, list[dict[str, Any]], dict[str, Any]]:
organ_df = pd.read_csv(organ_mapping_csv)
organ_required = {"token_id", "organ_id", "organ_label", "organ_weight"}
missing = sorted(organ_required - set(organ_df.columns))
if missing:
raise ValueError(f"{organ_mapping_csv} is missing required columns: {missing}")
organ_df = organ_df.copy()
organ_df["token_id"] = pd.to_numeric(organ_df["token_id"], errors="raise").astype(int)
organ_df["organ_weight"] = pd.to_numeric(
organ_df["organ_weight"], errors="raise"
).astype(float)
organ_df = organ_df[(organ_df["organ_id"].astype(str) != "") & (organ_df["organ_weight"] > 0)]
if organ_df.empty:
raise ValueError(f"{organ_mapping_csv} has no mapped organ rows.")
hfrs_df = pd.read_csv(hfrs_mapping_csv)
hfrs_required = {"token_id", "hfrs_weight"}
missing = sorted(hfrs_required - set(hfrs_df.columns))
if missing:
raise ValueError(f"{hfrs_mapping_csv} is missing required columns: {missing}")
hfrs_df = hfrs_df.copy()
hfrs_df["token_id"] = pd.to_numeric(hfrs_df["token_id"], errors="raise").astype(int)
hfrs_df["hfrs_weight"] = pd.to_numeric(
hfrs_df["hfrs_weight"], errors="raise"
).astype(float)
hfrs_df = hfrs_df[hfrs_df["hfrs_weight"] > 0]
if hfrs_df.empty:
raise ValueError(f"{hfrs_mapping_csv} has no non-zero HFRS weights.")
union_disease_ids = np.asarray(
sorted(
set(organ_df["token_id"].astype(int).tolist())
| set(hfrs_df["token_id"].astype(int).tolist())
),
dtype=np.int64,
)
union_pos = {int(token): i for i, token in enumerate(union_disease_ids.tolist())}
organ_ids = sorted(organ_df["organ_id"].astype(str).unique().tolist())
organ_pos = {organ_id: i for i, organ_id in enumerate(organ_ids)}
organ_matrix = np.zeros((len(organ_ids), union_disease_ids.size), dtype=np.float32)
organ_meta_by_id = {}
for _, row in organ_df.iterrows():
organ_id = str(row["organ_id"])
token = int(row["token_id"])
organ_matrix[organ_pos[organ_id], union_pos[token]] = 1.0
organ_meta_by_id.setdefault(
organ_id,
{
"index_type": "organ_involvement",
"index_id": organ_id,
"index_label": str(row["organ_label"]),
},
)
organ_meta = [organ_meta_by_id[organ_id] for organ_id in organ_ids]
hfrs_weights = np.zeros(union_disease_ids.size, dtype=np.float32)
for _, row in hfrs_df.iterrows():
hfrs_weights[union_pos[int(row["token_id"])]] = float(row["hfrs_weight"])
hfrs_meta = {
"index_type": "frailty_risk",
"index_id": "deephealth_hfrs",
"index_label": "DeepHealth-HFRS frailty risk index",
}
matrices = [
{
"kind": "organ_involvement",
"matrix": organ_matrix,
"meta": organ_meta,
},
{
"kind": "frailty_risk",
"weights": hfrs_weights,
"meta": hfrs_meta,
},
]
return union_disease_ids, matrices, {
"organ_mapped_tokens": int(organ_df["token_id"].nunique()),
"hfrs_mapped_tokens": int(hfrs_df["token_id"].nunique()),
}
def _config_split_indices(
n: int,
cfg: dict[str, Any],
eval_split: str,
subset_size: int,
) -> np.ndarray:
args = argparse.Namespace(
train_ratio=None,
val_ratio=None,
test_ratio=None,
seed=None,
eval_split=eval_split,
dataset_subset_size=subset_size if subset_size > 0 else None,
)
class _Sized:
def __len__(self) -> int:
return n
return make_eval_indices(_Sized(), args, cfg)
def _eligible_landmark_rows(
dataset: Any,
subset_indices: np.ndarray,
landmark_ages: np.ndarray,
*,
min_history_events: int,
) -> list[dict[str, Any]]:
rows: list[dict[str, Any]] = []
special = np.asarray([PAD_IDX, CHECKUP_IDX, NO_EVENT_IDX], dtype=np.int64)
for patient_id, dataset_index in enumerate(subset_indices.tolist()):
sample = dataset.samples[int(dataset_index)]
seq_event = np.asarray(sample["event_seq"], dtype=np.int64)
seq_time = np.asarray(sample["time_seq"], dtype=np.float32)
tgt_event = np.asarray(sample["target_event_seq"], dtype=np.int64)
tgt_time = np.asarray(sample["target_time_seq"], dtype=np.float32)
if seq_event.size == 0 or tgt_event.size == 0:
continue
full_event = np.concatenate([seq_event, tgt_event[-1:]])
full_time = np.concatenate([seq_time, tgt_time[-1:]])
followup_end = float(np.max(full_time))
for landmark_age in landmark_ages.tolist():
t_query = np.float32(float(landmark_age))
if not (followup_end > float(t_query)):
continue
prefix_mask = full_time <= t_query
if not np.any(prefix_mask):
continue
prefix_events = full_event[prefix_mask].astype(np.int64, copy=False)
valid_history = ~np.isin(prefix_events, special)
if int(valid_history.sum()) < int(min_history_events):
continue
rows.append(
{
"patient_id": int(patient_id),
"dataset_index": int(dataset_index),
"sex": int(sample["sex"]),
"landmark_age": t_query,
"t_query": t_query,
"followup_end_time": np.float32(followup_end),
"event_seq": prefix_events,
"time_seq": full_time[prefix_mask].astype(np.float32, copy=False),
"other_type": np.asarray(sample["other_type"], dtype=np.int64),
"other_value": np.asarray(sample["other_value"], dtype=np.float32),
"other_value_kind": np.asarray(sample["other_value_kind"], dtype=np.int64),
"other_time": np.asarray(sample["other_time"], dtype=np.float32),
}
)
return rows
def _row_to_worker_spec(row: dict[str, Any]) -> dict[str, Any]:
return {
"patient_id": int(row["patient_id"]),
"dataset_index": int(row["dataset_index"]),
"landmark_age": float(row["landmark_age"]),
"followup_end_time": float(row["followup_end_time"]),
}
def _materialize_worker_rows(
dataset: Any,
row_specs: list[dict[str, Any]],
) -> list[dict[str, Any]]:
rows: list[dict[str, Any]] = []
for spec in row_specs:
sample = dataset.samples[int(spec["dataset_index"])]
seq_event = np.asarray(sample["event_seq"], dtype=np.int64)
seq_time = np.asarray(sample["time_seq"], dtype=np.float32)
tgt_event = np.asarray(sample["target_event_seq"], dtype=np.int64)
tgt_time = np.asarray(sample["target_time_seq"], dtype=np.float32)
full_event = np.concatenate([seq_event, tgt_event[-1:]])
full_time = np.concatenate([seq_time, tgt_time[-1:]])
t_query = np.float32(float(spec["landmark_age"]))
prefix_mask = full_time <= t_query
rows.append(
{
"patient_id": int(spec["patient_id"]),
"dataset_index": int(spec["dataset_index"]),
"sex": int(sample["sex"]),
"landmark_age": t_query,
"t_query": t_query,
"followup_end_time": np.float32(float(spec["followup_end_time"])),
"event_seq": full_event[prefix_mask].astype(np.int64, copy=False),
"time_seq": full_time[prefix_mask].astype(np.float32, copy=False),
"other_type": np.asarray(sample["other_type"], dtype=np.int64),
"other_value": np.asarray(sample["other_value"], dtype=np.float32),
"other_value_kind": np.asarray(sample["other_value_kind"], dtype=np.int64),
"other_time": np.asarray(sample["other_time"], dtype=np.float32),
}
)
return rows
class HistoricalReadoutDataset(IterableDataset):
def __init__(self, rows: list[dict[str, Any]]) -> None:
super().__init__()
self.rows = rows
def __iter__(self) -> Iterable[dict[str, torch.Tensor]]:
worker = get_worker_info()
if worker is None:
start, step = 0, 1
else:
start, step = int(worker.id), int(worker.num_workers)
for row_idx in range(start, len(self.rows), step):
row = self.rows[row_idx]
grid = build_readout_grid(
event_seq=row["event_seq"],
time_seq=row["time_seq"],
other_type=row["other_type"],
other_time=row["other_time"],
t_query=float(row["t_query"]),
)
if grid.size == 0:
continue
end_times = np.concatenate([grid[1:], np.asarray([row["t_query"]], dtype=np.float32)])
deltas = np.maximum(end_times - grid, 0.0).astype(np.float32)
valid = deltas > 0
for query_time, delta in zip(grid[valid].tolist(), deltas[valid].tolist()):
yield _make_readout_job(row, row_idx, query_time, delta)
def _make_readout_job(
row: dict[str, Any],
row_idx: int,
query_time: float,
delta: float,
) -> dict[str, torch.Tensor]:
return {
"event_seq": torch.from_numpy(np.asarray(row["event_seq"], dtype=np.int64)).long(),
"time_seq": torch.from_numpy(np.asarray(row["time_seq"], dtype=np.float32)).float(),
"sex": torch.tensor(int(row["sex"]), dtype=torch.long),
"other_type": torch.from_numpy(np.asarray(row["other_type"], dtype=np.int64)).long(),
"other_value": torch.from_numpy(np.asarray(row["other_value"], dtype=np.float32)).float(),
"other_value_kind": torch.from_numpy(
np.asarray(row["other_value_kind"], dtype=np.int64)
).long(),
"other_time": torch.from_numpy(np.asarray(row["other_time"], dtype=np.float32)).float(),
"query_time": torch.tensor(float(query_time), dtype=torch.float32),
"delta": torch.tensor(float(delta), dtype=torch.float32),
"row_idx": torch.tensor(int(row_idx), dtype=torch.long),
}
def _collate_readout_jobs(batch: list[dict[str, torch.Tensor]]) -> dict[str, torch.Tensor]:
event_seq = pad_sequence(
[x["event_seq"] for x in batch], batch_first=True, padding_value=PAD_IDX
)
return {
"event_seq": event_seq,
"time_seq": pad_sequence(
[x["time_seq"] for x in batch], batch_first=True, padding_value=0.0
),
"padding_mask": event_seq > PAD_IDX,
"sex": torch.stack([x["sex"] for x in batch]),
"other_type": pad_sequence(
[x["other_type"] for x in batch], batch_first=True, padding_value=0
),
"other_value": pad_sequence(
[x["other_value"] for x in batch], batch_first=True, padding_value=0.0
),
"other_value_kind": pad_sequence(
[x["other_value_kind"] for x in batch], batch_first=True, padding_value=0
),
"other_time": pad_sequence(
[x["other_time"] for x in batch], batch_first=True, padding_value=0.0
),
"query_time": torch.stack([x["query_time"] for x in batch]),
"delta": torch.stack([x["delta"] for x in batch]),
"row_idx": torch.stack([x["row_idx"] for x in batch]),
}
@torch.inference_mode()
def _readout_probabilities(
*,
ctx: Any,
batch: dict[str, torch.Tensor],
disease_ids: np.ndarray,
) -> torch.Tensor:
event = batch["event_seq"].long().to(ctx.device, non_blocking=True)
hidden = ctx.model(
event_seq=event,
time_seq=batch["time_seq"].float().to(ctx.device, non_blocking=True),
sex=batch["sex"].long().to(ctx.device, non_blocking=True),
padding_mask=event > PAD_IDX,
t_query=batch["query_time"].float().to(ctx.device, non_blocking=True),
other_type=batch["other_type"].long().to(ctx.device, non_blocking=True),
other_value=batch["other_value"].float().to(ctx.device, non_blocking=True),
other_value_kind=batch["other_value_kind"].long().to(ctx.device, non_blocking=True),
other_time=batch["other_time"].float().to(ctx.device, non_blocking=True),
target_mode="all_future",
)
deltas = batch["delta"].detach().cpu().numpy().astype(np.float32, copy=False)
prob = probabilities_from_hidden(
ctx=ctx,
hidden=hidden,
disease_ids=disease_ids,
deltas=deltas,
)
return torch.as_tensor(prob, dtype=torch.float32, device=ctx.device)
def _project_rows(
*,
rows: list[dict[str, Any]],
survival_by_row: torch.Tensor,
matrices: list[dict[str, Any]],
ctx: Any,
) -> list[dict[str, Any]]:
disease_expression = 1.0 - survival_by_row.clamp(0.0, 1.0)
disease_intensity = -torch.log(survival_by_row.clamp(1e-7, 1.0))
out: list[dict[str, Any]] = []
organ_matrix = torch.as_tensor(
matrices[0]["matrix"], dtype=torch.float32, device=ctx.device
)
organ_values = -torch.expm1(-torch.matmul(disease_intensity, organ_matrix.T))
organ_values_np = organ_values.detach().cpu().numpy()
hfrs_weights = torch.as_tensor(
matrices[1]["weights"], dtype=torch.float32, device=ctx.device
)
hfrs_values = torch.matmul(disease_expression, hfrs_weights)
hfrs_values_np = hfrs_values.detach().cpu().numpy()
for row_idx, row in enumerate(rows):
base = {
"patient_id": row["patient_id"],
"dataset_index": row["dataset_index"],
"sex": row["sex"],
"landmark_age": float(row["landmark_age"]),
"t_query": float(row["t_query"]),
"followup_end_time": float(row["followup_end_time"]),
}
for dim_idx, meta in enumerate(matrices[0]["meta"]):
out.append(
{
**base,
"index_type": meta["index_type"],
"index_id": meta["index_id"],
"index_label": meta["index_label"],
"index_value": float(organ_values_np[row_idx, dim_idx]),
}
)
out.append(
{
**base,
"index_type": matrices[1]["meta"]["index_type"],
"index_id": matrices[1]["meta"]["index_id"],
"index_label": matrices[1]["meta"]["index_label"],
"index_value": float(hfrs_values_np[row_idx]),
}
)
return out
def _compute_rows(
*,
rows: list[dict[str, Any]],
disease_ids: np.ndarray,
matrices: list[dict[str, Any]],
readout_batch_size: int,
num_workers: int,
ctx: Any,
log_prefix: str,
) -> tuple[list[dict[str, Any]], int, dict[str, float]]:
survival_by_row = torch.ones(
(len(rows), disease_ids.size), dtype=torch.float32, device=ctx.device
)
loader = DataLoader(
HistoricalReadoutDataset(rows),
batch_size=max(1, int(readout_batch_size)),
collate_fn=_collate_readout_jobs,
num_workers=max(0, int(num_workers)),
pin_memory=ctx.device.type == "cuda",
persistent_workers=int(num_workers) > 0,
prefetch_factor=2 if int(num_workers) > 0 else None,
)
readout_jobs = 0
n_batches = 0
forward_sec = 0.0
reduce_sec = 0.0
for batch in loader:
t0 = time.perf_counter()
prob = _readout_probabilities(ctx=ctx, batch=batch, disease_ids=disease_ids)
if ctx.device.type == "cuda":
torch.cuda.synchronize(ctx.device)
forward_sec += time.perf_counter() - t0
t1 = time.perf_counter()
row_indices = batch["row_idx"].long().to(ctx.device, non_blocking=True)
interval_survival = 1.0 - prob.clamp(0.0, 1.0)
if hasattr(survival_by_row, "scatter_reduce_"):
survival_by_row.scatter_reduce_(
dim=0,
index=row_indices[:, None].expand_as(interval_survival),
src=interval_survival,
reduce="prod",
include_self=True,
)
else:
for job_idx in range(interval_survival.shape[0]):
survival_by_row[int(row_indices[job_idx].item())] *= interval_survival[job_idx]
if ctx.device.type == "cuda":
torch.cuda.synchronize(ctx.device)
reduce_sec += time.perf_counter() - t1
n_batches += 1
readout_jobs += int(batch["row_idx"].numel())
if n_batches == 1 or n_batches % 50 == 0:
print(
f"{log_prefix} processed {readout_jobs} readout jobs in {n_batches} batches",
flush=True,
)
t2 = time.perf_counter()
out = _project_rows(
rows=rows,
survival_by_row=survival_by_row,
matrices=matrices,
ctx=ctx,
)
if ctx.device.type == "cuda":
torch.cuda.synchronize(ctx.device)
reduce_sec += time.perf_counter() - t2
return out, readout_jobs, {"forward_sec": forward_sec, "reduce_sec": reduce_sec}
def _write_rows_csv(rows: list[dict[str, Any]], output_path: Path) -> int:
df = pd.DataFrame(rows)
df.to_csv(output_path, index=False)
return int(len(df))
def _concat_csv_shards(shard_paths: list[Path], output_path: Path) -> None:
wrote_header = False
with output_path.open("w", encoding="utf-8", newline="") as out_f:
for shard_path in shard_paths:
with shard_path.open("r", encoding="utf-8", newline="") as in_f:
header = in_f.readline()
if not wrote_header:
out_f.write(header)
wrote_header = True
for line in in_f:
out_f.write(line)
shard_path.unlink(missing_ok=True)
def _estimate_jobs(row: dict[str, Any]) -> int:
grid = build_readout_grid(
event_seq=row["event_seq"],
time_seq=row["time_seq"],
other_type=row["other_type"],
other_time=row["other_time"],
t_query=float(row["t_query"]),
)
if grid.size == 0:
return 1
end_times = np.concatenate([grid[1:], np.asarray([row["t_query"]], dtype=np.float32)])
return max(int(np.sum(np.maximum(end_times - grid, 0.0) > 0)), 1)
def _split_rows(rows: list[dict[str, Any]], devices: list[str | None]) -> list[tuple[str | None, list[dict[str, Any]]]]:
if len(devices) <= 1:
return [(devices[0], rows)]
buckets: list[list[dict[str, Any]]] = [[] for _ in devices]
loads = np.zeros(len(devices), dtype=np.int64)
for row in sorted(rows, key=_estimate_jobs, reverse=True):
idx = int(np.argmin(loads))
buckets[idx].append(row)
loads[idx] += _estimate_jobs(row)
return [(device, bucket) for device, bucket in zip(devices, buckets) if bucket]
def _worker(payload: dict[str, Any]) -> dict[str, Any]:
device = payload["device"]
shard_path = Path(payload["shard_path"])
print(f"[Index worker {device}] starting with {len(payload['row_specs'])} rows", flush=True)
ctx = load_deephealth_context(payload["run_path"], device=device)
rows = _materialize_worker_rows(ctx.dataset, payload["row_specs"])
out, readout_jobs, timings = _compute_rows(
rows=rows,
disease_ids=payload["disease_ids"],
matrices=payload["matrices"],
readout_batch_size=int(payload["readout_batch_size"]),
num_workers=int(payload["num_workers"]),
ctx=ctx,
log_prefix=f"[Index worker {device}]",
)
t0 = time.perf_counter()
row_count = _write_rows_csv(out, shard_path)
timings["write_csv_sec"] = time.perf_counter() - t0
print(f"[Index worker {device}] wrote {row_count} rows to {shard_path}", flush=True)
return {
"shard_path": str(shard_path),
"row_count": row_count,
"readout_jobs": readout_jobs,
"timings": timings,
}
def main() -> None:
parser = argparse.ArgumentParser(
description="Compute DeepHealth organ involvement and frailty risk indices."
)
parser.add_argument("--run_path", type=str, required=True)
parser.add_argument(
"--organ_mapping_csv",
type=str,
default="organ_involvement_label_mapping.csv",
)
parser.add_argument("--hfrs_mapping_csv", type=str, default="uk_hfrs_label_mapping.csv")
parser.add_argument("--output_path", type=str, default=None)
parser.add_argument(
"--eval_split",
type=str,
default="test",
choices=["train", "val", "valid", "validation", "test", "all"],
)
parser.add_argument("--landmark_ages", type=str, default=None)
parser.add_argument("--landmark_start", type=float, default=40.0)
parser.add_argument("--landmark_stop", type=float, default=80.0)
parser.add_argument("--landmark_step", type=float, default=5.0)
parser.add_argument("--min_history_events", type=int, default=1)
parser.add_argument("--dataset_subset_size", type=int, default=0)
parser.add_argument("--readout_batch_size", type=int, default=8192)
parser.add_argument("--num_workers", type=int, default=4)
parser.add_argument("--device", type=str, default=None)
parser.add_argument("--devices", type=str, default=None)
parser.add_argument(
"--mp_start_method",
type=str,
default="fork",
choices=["fork", "forkserver"],
)
args = parser.parse_args()
run_path = Path(args.run_path)
devices = _parse_devices(args)
initial_device = "cpu" if len(devices) > 1 else devices[0]
ctx = load_deephealth_context(run_path, device=initial_device)
disease_ids, matrices, mapping_counts = _load_index_matrices(
organ_mapping_csv=Path(args.organ_mapping_csv),
hfrs_mapping_csv=Path(args.hfrs_mapping_csv),
)
landmark_ages = _parse_landmark_ages(args)
eval_split = str(args.eval_split).lower()
if eval_split in {"valid", "validation"}:
eval_split = "val"
subset_indices = _config_split_indices(
len(ctx.dataset),
ctx.cfg,
eval_split,
int(args.dataset_subset_size),
)
rows = _eligible_landmark_rows(
ctx.dataset,
subset_indices,
landmark_ages,
min_history_events=int(args.min_history_events),
)
if not rows:
raise RuntimeError("No eligible landmark rows.")
output_path = Path(args.output_path) if args.output_path else (
run_path / f"deephealth_indices_{eval_split}.csv"
)
output_path.parent.mkdir(parents=True, exist_ok=True)
chunks = _split_rows(rows, devices)
for device, chunk in chunks:
print(
f"Assigned {len(chunk)} rows / ~{sum(_estimate_jobs(r) for r in chunk)} "
f"readout jobs to {device}",
flush=True,
)
total_readout_jobs = 0
timings = {"forward_sec": 0.0, "reduce_sec": 0.0, "write_csv_sec": 0.0}
if len(chunks) == 1:
out, total_readout_jobs, chunk_timings = _compute_rows(
rows=rows,
disease_ids=disease_ids,
matrices=matrices,
readout_batch_size=int(args.readout_batch_size),
num_workers=int(args.num_workers),
ctx=ctx,
log_prefix="[Index main]",
)
for key, value in chunk_timings.items():
timings[key] += float(value)
t0 = time.perf_counter()
_write_rows_csv(out, output_path)
timings["write_csv_sec"] = time.perf_counter() - t0
else:
del ctx
payloads = [
{
"device": device,
"run_path": str(run_path),
"shard_path": str(
output_path.with_name(
f"{output_path.stem}.part{part_idx:03d}{output_path.suffix}"
)
),
"row_specs": [_row_to_worker_spec(row) for row in chunk],
"disease_ids": disease_ids,
"matrices": matrices,
"readout_batch_size": int(args.readout_batch_size),
"num_workers": int(args.num_workers),
}
for part_idx, (device, chunk) in enumerate(chunks)
]
shard_paths = []
with ProcessPoolExecutor(
max_workers=len(payloads),
mp_context=mp.get_context(args.mp_start_method),
) as executor:
futures = [executor.submit(_worker, payload) for payload in payloads]
for future in tqdm(
as_completed(futures),
total=len(futures),
desc="Computing DeepHealth index chunks",
dynamic_ncols=True,
):
result = future.result()
shard_paths.append(Path(result["shard_path"]))
total_readout_jobs += int(result["readout_jobs"])
for key, value in result["timings"].items():
timings[key] += float(value)
t0 = time.perf_counter()
_concat_csv_shards(sorted(shard_paths), output_path)
timings["write_csv_sec"] += time.perf_counter() - t0
print(f"Run path: {run_path}")
print(f"Eval split: {eval_split}")
print(f"Landmark rows: {len(rows)}")
print(f"Readout jobs: {total_readout_jobs}")
print(f"Union disease tokens: {disease_ids.size}")
print(f"Organ mapped tokens: {mapping_counts['organ_mapped_tokens']}")
print(f"HFRS mapped tokens: {mapping_counts['hfrs_mapped_tokens']}")
print(
"Timing seconds: "
f"forward={timings['forward_sec']:.2f}, "
f"reduce={timings['reduce_sec']:.2f}, "
f"write_csv={timings['write_csv_sec']:.2f}"
)
print(f"Devices: {', '.join(str(d) for d, _ in chunks)}")
print(f"Output: {output_path}")
if __name__ == "__main__":
main()

File diff suppressed because it is too large Load Diff

View File

@@ -1,214 +0,0 @@
from __future__ import annotations
import argparse
from pathlib import Path
try:
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
except ModuleNotFoundError as exc:
raise ModuleNotFoundError(
"plot_deephealth_index_trajectories.py requires matplotlib. "
"Install it in the server environment before running this script."
) from exc
import numpy as np
import pandas as pd
REQUIRED_COLUMNS = {
"patient_id",
"sex",
"landmark_age",
"index_type",
"index_id",
"index_label",
"index_value",
}
def _sex_label(value: object) -> str:
text = str(value).strip().lower()
if text in {"0", "0.0", "female", "f", "woman"}:
return "female"
if text in {"1", "1.0", "male", "m", "man"}:
return "male"
return text or "unknown"
def _load_index_csv(path: Path) -> pd.DataFrame:
header = pd.read_csv(path, nrows=0)
missing = sorted(REQUIRED_COLUMNS - set(header.columns))
if missing:
raise ValueError(f"{path} is missing required columns: {missing}")
df = pd.read_csv(path, usecols=sorted(REQUIRED_COLUMNS))
df["sex_label"] = df["sex"].map(_sex_label)
df["landmark_age"] = pd.to_numeric(df["landmark_age"], errors="raise")
df["index_value"] = pd.to_numeric(df["index_value"], errors="raise")
return df
def _sample_patients(df: pd.DataFrame, *, n_per_sex: int, seed: int) -> dict[str, np.ndarray]:
rng = np.random.default_rng(int(seed))
samples: dict[str, np.ndarray] = {}
for sex_label in ["female", "male"]:
ids = np.asarray(sorted(df.loc[df["sex_label"] == sex_label, "patient_id"].unique()))
if ids.size == 0:
samples[sex_label] = np.asarray([], dtype=np.int64)
continue
samples[sex_label] = np.asarray(
rng.choice(ids, size=min(int(n_per_sex), int(ids.size)), replace=False)
)
return samples
def _plot_sampled_trajectories(
df: pd.DataFrame,
*,
index_type: str,
selected: dict[str, np.ndarray],
output_dir: Path,
) -> None:
sub = df[df["index_type"] == index_type].copy()
if sub.empty:
return
if index_type == "organ_involvement":
total = (
sub.groupby(["patient_id", "sex_label", "landmark_age"], as_index=False)[
"index_value"
]
.mean()
.rename(columns={"index_value": "trajectory_value"})
)
title = "mean organ involvement"
else:
total = sub.rename(columns={"index_value": "trajectory_value"})
title = "frailty risk"
mean_df = (
total.groupby(["sex_label", "landmark_age"], as_index=False)["trajectory_value"]
.mean()
.sort_values("landmark_age")
)
fig, ax = plt.subplots(figsize=(9.5, 5.5), dpi=160)
colors = {"female": "#b83280", "male": "#2563eb"}
for sex_label in ["female", "male"]:
for pid in selected.get(sex_label, np.asarray([], dtype=np.int64)):
one = total[total["patient_id"] == pid].sort_values("landmark_age")
if one.empty:
continue
ax.plot(
one["landmark_age"],
one["trajectory_value"],
color=colors.get(sex_label, "0.4"),
alpha=0.22,
linewidth=1.2,
)
mean_one = mean_df[mean_df["sex_label"] == sex_label]
if not mean_one.empty:
ax.plot(
mean_one["landmark_age"],
mean_one["trajectory_value"],
color=colors.get(sex_label, "0.4"),
linewidth=2.6,
label=f"{sex_label} mean",
)
ax.set_title(f"{index_type}: sampled trajectories and sex-specific means")
ax.set_xlabel("Landmark age")
ax.set_ylabel(title)
ax.grid(True, alpha=0.25)
ax.legend(frameon=False)
fig.tight_layout()
fig.savefig(output_dir / f"{index_type}_sampled_trajectories_by_sex.png")
plt.close(fig)
def _plot_top_dimensions(
df: pd.DataFrame,
*,
output_dir: Path,
top_n: int,
) -> None:
sub = df[df["index_type"] == "organ_involvement"].copy()
if sub.empty:
return
order = (
sub.groupby("index_id")["index_value"]
.mean()
.sort_values(ascending=False)
.head(int(top_n))
.index.tolist()
)
n = len(order)
if n == 0:
return
ncols = min(3, n)
nrows = int(np.ceil(n / ncols))
fig, axes = plt.subplots(nrows, ncols, figsize=(4.0 * ncols, 3.0 * nrows), dpi=160)
axes_arr = np.asarray(axes).reshape(-1)
colors = {"female": "#b83280", "male": "#2563eb"}
for ax, index_id in zip(axes_arr, order):
one = sub[sub["index_id"] == index_id]
mean_df = (
one.groupby(["sex_label", "landmark_age"], as_index=False)["index_value"]
.mean()
.sort_values("landmark_age")
)
for sex_label in ["female", "male"]:
m = mean_df[mean_df["sex_label"] == sex_label]
if not m.empty:
ax.plot(
m["landmark_age"],
m["index_value"],
color=colors.get(sex_label, "0.4"),
linewidth=1.8,
label=sex_label,
)
ax.set_title(str(index_id), fontsize=9)
ax.set_xlabel("Age")
ax.set_ylabel("Organ involvement")
ax.grid(True, alpha=0.22)
for ax in axes_arr[n:]:
ax.axis("off")
handles, labels = axes_arr[0].get_legend_handles_labels()
if handles:
fig.legend(handles, labels, loc="upper right", frameon=False)
fig.suptitle("organ_involvement: sex-specific mean trajectories for top dimensions")
fig.tight_layout(rect=(0, 0, 0.98, 0.96))
fig.savefig(output_dir / "organ_involvement_top_dimensions_by_sex.png")
plt.close(fig)
def main() -> None:
parser = argparse.ArgumentParser(
description="Plot DeepHealth organ involvement and frailty risk trajectories."
)
parser.add_argument("--input_csv", type=str, required=True)
parser.add_argument("--output_dir", type=str, default=None)
parser.add_argument("--n_per_sex", type=int, default=10)
parser.add_argument("--seed", type=int, default=2026)
parser.add_argument("--top_n", type=int, default=12)
args = parser.parse_args()
input_csv = Path(args.input_csv)
output_dir = Path(args.output_dir) if args.output_dir else input_csv.parent
output_dir.mkdir(parents=True, exist_ok=True)
df = _load_index_csv(input_csv)
selected = _sample_patients(df, n_per_sex=int(args.n_per_sex), seed=int(args.seed))
for index_type in ["organ_involvement", "frailty_risk"]:
_plot_sampled_trajectories(
df,
index_type=index_type,
selected=selected,
output_dir=output_dir,
)
_plot_top_dimensions(df, output_dir=output_dir, top_n=int(args.top_n))
print(f"Input: {input_csv}")
print(f"Output directory: {output_dir}")
if __name__ == "__main__":
main()

File diff suppressed because it is too large Load Diff

View File

@@ -1,53 +0,0 @@
hfrs_source_code,hfrs_weight,missing_reason,hfrs_source
R00,0.7,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R02,1.0,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R11,0.3,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R13,0.8,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R26,2.6,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R29,3.6,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R31,3.0,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R32,1.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R33,1.3,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R40,2.5,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R41,2.7,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R44,1.6,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R45,1.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R47,1.0,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R50,0.1,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R54,2.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R55,1.8,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R56,2.6,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R63,0.9,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R69,1.3,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R79,0.6,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
R94,1.4,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S00,3.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S01,1.1,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S06,2.4,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S09,1.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S22,1.8,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S32,1.4,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S42,2.3,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S51,0.5,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S72,1.4,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
S80,2.0,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
T83,2.4,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
U80,0.8,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
W01,0.9,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
W06,1.1,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
W10,0.9,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
W18,2.1,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
W19,3.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
X59,1.5,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Y84,0.7,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Y95,1.2,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z22,1.7,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z50,2.1,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z60,1.8,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z73,0.6,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z74,1.1,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z75,2.0,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z87,1.5,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z91,0.5,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z93,1.0,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
Z99,0.8,not_present_in_labels_csv,Gilbert et al. Lancet 2018 supplementary appendix Table A2
1 hfrs_source_code hfrs_weight missing_reason hfrs_source
2 R00 0.7 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
3 R02 1.0 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
4 R11 0.3 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
5 R13 0.8 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
6 R26 2.6 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
7 R29 3.6 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
8 R31 3.0 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
9 R32 1.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
10 R33 1.3 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
11 R40 2.5 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
12 R41 2.7 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
13 R44 1.6 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
14 R45 1.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
15 R47 1.0 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
16 R50 0.1 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
17 R54 2.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
18 R55 1.8 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
19 R56 2.6 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
20 R63 0.9 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
21 R69 1.3 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
22 R79 0.6 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
23 R94 1.4 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
24 S00 3.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
25 S01 1.1 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
26 S06 2.4 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
27 S09 1.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
28 S22 1.8 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
29 S32 1.4 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
30 S42 2.3 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
31 S51 0.5 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
32 S72 1.4 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
33 S80 2.0 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
34 T83 2.4 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
35 U80 0.8 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
36 W01 0.9 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
37 W06 1.1 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
38 W10 0.9 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
39 W18 2.1 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
40 W19 3.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
41 X59 1.5 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
42 Y84 0.7 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
43 Y95 1.2 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
44 Z22 1.7 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
45 Z50 2.1 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
46 Z60 1.8 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
47 Z73 0.6 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
48 Z74 1.1 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
49 Z75 2.0 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
50 Z87 1.5 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
51 Z91 0.5 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
52 Z93 1.0 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2
53 Z99 0.8 not_present_in_labels_csv Gilbert et al. Lancet 2018 supplementary appendix Table A2