diff --git a/.gitignore b/.gitignore index 07523f7..c4c85e6 100644 --- a/.gitignore +++ b/.gitignore @@ -89,6 +89,13 @@ dump.rdb "-w" "]" +# 临时/调试脚本与一次性验证输出(不进版本库) +backend/scripts/_tmp_*.py +backend/scripts/_apply_trait_seed_sql.py +backend/scripts/breeding_stats/_debug_*.py +backend/scripts/breeding_stats/_sparse_golden_check.py +backend/scripts/breeding_stats/_golden_out.txt + # 本地二进制/运行环境(不提交) redis/ backend/.venv.broken/ \ No newline at end of file diff --git a/backend/app/api/v1/module_bre/statistics/service.py b/backend/app/api/v1/module_bre/statistics/service.py index 0131528..88f7de0 100644 --- a/backend/app/api/v1/module_bre/statistics/service.py +++ b/backend/app/api/v1/module_bre/statistics/service.py @@ -26,7 +26,8 @@ import numpy as np from datetime import date, datetime from typing import Any -from sqlalchemy import String, case, func, or_, select, text, update +from sqlalchemy import Integer, String, case, cast, false, func, or_, select, text, update +from sqlalchemy.dialects.postgresql import ARRAY from sqlalchemy.ext.asyncio import AsyncSession from app.core.base_schema import AuthSchema @@ -84,6 +85,19 @@ def _now() -> datetime: return datetime.now() +def _id_in(column, ids): + """按 IN (SELECT unnest(:arr)) 绑定大 id 列表为单个数组参数。 + + asyncpg 单条查询总绑定参数上限 32767,常规 column.in_(ids) 会在 + bre_tree 等表 id 规模过万时直接崩(IN 子句分块无济于事,因为是全查询 + 计数)。空列表返回恒假条件,避免 SQLAlchemy 渲染出畸形的 IN (NULL)。 + """ + ids = list(ids) + if not ids: + return false() + return column.in_(select(func.unnest(cast(ids, ARRAY(Integer))))) + + # 交配设计(配合力 GCA/SCA 的算法与 ANOVA 自由度随设计而异;与 combining.solve 分支一致) _COMBINING_DESIGNS = combining.DESIGNS _DESIGN_LABELS = { @@ -438,6 +452,27 @@ class StatisticsService: )).all() return {int(r.tree_id) for r in rows} + @staticmethod + def _germ_ancestor_closure(tree_germs: set[int], + g_parent_of: dict[int, tuple[int | None, int | None]]) -> set[int]: + """种质级祖先闭包:从选中树直系亲本向上取可达祖先(bre_pedigree 只追祖先链)。 + + 旧实现把全库 bre_pedigree 子代种质无条件并进模型(选中 400 树 → n≈1204),无关种质 + 在 MME 中零贡献(无表型、与选中树无亲缘),EBV/h² 逐位不变,但 n=树数+~800 一越过 + N_SPARSE 就误入稀疏精确迹二分(~11min/次)。收口到真实系谱规模后,常见「选几百棵 + 候选树」落稠密秒级。与 _germplasm_pedigree 的多源 BFS 思路一致。 + """ + out = set(tree_germs) + stack = list(tree_germs) + while stack: + g = stack.pop() + gd, gs = g_parent_of.get(g, (None, None)) + for pp in (gd, gs): + if pp is not None and pp not in out: + out.add(pp) + stack.append(pp) + return out + async def _build_pedigree(self, tree_ids: list[int]) -> list[dict]: """为给定树集合构建系谱(树=子代、种质=亲本,t/g 前缀;含种质级祖先 bre_pedigree)。 @@ -456,7 +491,7 @@ class StatisticsService: CrossCombinationModel.male_parent_id, ) .outerjoin(CrossCombinationModel, TreeModel.combination_id == CrossCombinationModel.id) - .where(TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + .where(_id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() pedigree: list[dict] = [] tree_germs: set[int] = set() @@ -489,13 +524,7 @@ class StatisticsService: cid = code2g.get(child_code) if cid is not None and (gd is not None or gs is not None): g_parent_of[cid] = (gd, gs) - all_germ_ids = set(tree_germs) - for g, (gd, gs) in g_parent_of.items(): - all_germ_ids.add(g) - if gd is not None: - all_germ_ids.add(gd) - if gs is not None: - all_germ_ids.add(gs) + all_germ_ids = self._germ_ancestor_closure(tree_germs, g_parent_of) for g in sorted(all_germ_ids): if g in g_parent_of: gd, gs = g_parent_of[g] @@ -686,7 +715,7 @@ class StatisticsService: TreeModel.germplasm_id, ) .outerjoin(CrossCombinationModel, TreeModel.combination_id == CrossCombinationModel.id) - .where(TreeModel.id.in_(tree_ids)) + .where(_id_in(TreeModel.id, tree_ids)) )).all() pedigree: list[dict] = [] tree_germs: set[int] = set() @@ -739,7 +768,7 @@ class StatisticsService: pos_rows = (await self.db.execute( select(TreeModel.id, TreeModel.plot_id, TreeModel.block_no, TreeModel.row_no, TreeModel.col_no) - .where(TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + .where(_id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() grids: dict[tuple[int, int], list[tuple[int, int, int]]] = {} for r in pos_rows: @@ -776,13 +805,7 @@ class StatisticsService: cid = code2g.get(child_code) if cid is not None and (gd is not None or gs is not None): g_parent_of[cid] = (gd, gs) - all_germ_ids = set(tree_germs) - for g, (gd, gs) in g_parent_of.items(): - all_germ_ids.add(g) - if gd is not None: - all_germ_ids.add(gd) - if gs is not None: - all_germ_ids.add(gs) + all_germ_ids = self._germ_ancestor_closure(tree_germs, g_parent_of) for g in sorted(all_germ_ids): if g in g_parent_of: gd, gs = g_parent_of[g] @@ -846,7 +869,7 @@ class StatisticsService: TreeModel.is_deleted.is_(False), # 与 site 分支对齐:只取本批次选中树,避免他树混入后 # clone_of/combo_of 缺失 → 归入虚构基因型组 f{None} - TreeModel.id.in_(tree_ids), + _id_in(TreeModel.id, tree_ids), ) .group_by(yr_obs.tree_id, yr_obs.evaluate_year) ) @@ -1093,7 +1116,7 @@ class StatisticsService: phenos_keyed = {f"t{tid}": v for tid, v in phenotypes.items()} bpos_rows = (await self.db.execute( select(TreeModel.id, TreeModel.block_no) - .where(TreeModel.id.in_(tree_ids), TreeModel.block_no.isnot(None), + .where(_id_in(TreeModel.id, tree_ids), TreeModel.block_no.isnot(None), TreeModel.is_deleted.is_(False)) )).all() block_of: dict[str, int] = {} @@ -1171,7 +1194,7 @@ class StatisticsService: if spatial: pos_rows = (await self.db.execute( select(TreeModel.id, TreeModel.row_no, TreeModel.col_no) - .where(TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + .where(_id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() coords: dict[str, tuple[int, int]] = {} missing_coord_ids: list[int] = [] @@ -1391,7 +1414,7 @@ class StatisticsService: TreeModel.germplasm_id, ) .outerjoin(CrossCombinationModel, TreeModel.combination_id == CrossCombinationModel.id) - .where(TreeModel.id.in_(tree_ids)) + .where(_id_in(TreeModel.id, tree_ids)) )).all() pedigree: list[dict] = [] tree_germs: set[int] = set() @@ -1434,7 +1457,7 @@ class StatisticsService: pos_rows = (await self.db.execute( select(TreeModel.id, TreeModel.plot_id, TreeModel.block_no, TreeModel.row_no, TreeModel.col_no) - .where(TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + .where(_id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() grids: dict[tuple[int, int], list[tuple[int, int, int]]] = {} for r in pos_rows: @@ -1470,13 +1493,7 @@ class StatisticsService: cid = code2g.get(child_code) if cid is not None and (gd is not None or gs is not None): g_parent_of[cid] = (gd, gs) - all_germ_ids = set(tree_germs) - for g, (gd, gs) in g_parent_of.items(): - all_germ_ids.add(g) - if gd is not None: - all_germ_ids.add(gd) - if gs is not None: - all_germ_ids.add(gs) + all_germ_ids = self._germ_ancestor_closure(tree_germs, g_parent_of) for g in sorted(all_germ_ids): if g in g_parent_of: gd, gs = g_parent_of[g] @@ -1766,7 +1783,7 @@ class StatisticsService: rel = {r.tree_id: float(r.reliability or 0.0) for r in rows} trees = (await self.db.execute( select(TreeModel.id, TreeModel.trial_study_id, TreeModel.combination_id) - .where(TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + .where(_id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() site_of = {t.id: t.trial_study_id for t in trees} combo_of = {t.id: t.combination_id for t in trees} @@ -2015,7 +2032,7 @@ class StatisticsService: if tree_ids: rows = (await self.db.execute( select(TreeModel.id, TreeModel.clone_id).where( - TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False) + _id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False) ) )).all() tree_map = {t.id: t.clone_id for t in rows} @@ -2023,7 +2040,7 @@ class StatisticsService: clone_map: dict[int, CloneModel] = {} if clone_ids: clones = (await self.db.execute( - select(CloneModel).where(CloneModel.id.in_(clone_ids), CloneModel.is_deleted.is_(False)) + select(CloneModel).where(_id_in(CloneModel.id, clone_ids), CloneModel.is_deleted.is_(False)) )).scalars().all() clone_map = {c.id: c for c in clones} agg: dict[int, dict[str, float | int]] = {} @@ -2533,7 +2550,7 @@ class StatisticsService: pid2code = {pid: c for c, pid in batch_ids.items()} pv_rows = (await self.db.execute( select(PredictionValueModel.prediction_id, PredictionValueModel.tree_id, PredictionValueModel.reliability).where( - PredictionValueModel.prediction_id.in_(list(batch_ids.values())), + _id_in(PredictionValueModel.prediction_id, list(batch_ids.values())), PredictionValueModel.tree_id.isnot(None), ) )).all() @@ -2677,7 +2694,7 @@ class StatisticsService: """ rows = (await self.db.execute( select(TreeModel.id, TreeModel.clone_id, TreeModel.tree_no).where( - TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False) + _id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False) ) )).all() clone_of = {t.id: t.clone_id for t in rows} @@ -2687,7 +2704,7 @@ class StatisticsService: if clone_ids: clones = (await self.db.execute( select(CloneModel.id, CloneModel.clone_code).where( - CloneModel.id.in_(clone_ids), CloneModel.is_deleted.is_(False) + _id_in(CloneModel.id, clone_ids), CloneModel.is_deleted.is_(False) ) )).all() clone_label = {c.id: c.clone_code for c in clones} @@ -3011,7 +3028,7 @@ class StatisticsService: if refs: pv_rows = (await self.db.execute( select(PredictionValueModel.tree_id, PredictionValueModel.reliability).where( - PredictionValueModel.prediction_id.in_(list(refs.values())), + _id_in(PredictionValueModel.prediction_id, list(refs.values())), PredictionValueModel.tree_id.isnot(None), ) )).all() @@ -3048,11 +3065,11 @@ class StatisticsService: all_tids = {t for u in units for t in u["tree_ids"]} tree_rows = (await self.db.execute( select(TreeModel.id, TreeModel.combination_id, TreeModel.tree_no, TreeModel.stage) - .where(TreeModel.id.in_(all_tids)) + .where(_id_in(TreeModel.id, all_tids)) )).all() tree_map = {t.id: t for t in tree_rows} existing = (await self.db.execute( - select(SelectionResultModel).where(SelectionResultModel.tree_id.in_(all_tids)) + select(SelectionResultModel).where(_id_in(SelectionResultModel.tree_id, all_tids)) )).scalars().all() existing_by_tree = {e.tree_id: e for e in existing} @@ -3144,7 +3161,7 @@ class StatisticsService: """ stmt = select(TreeModel.id, TreeModel.dam_id, TreeModel.sire_id).where(TreeModel.is_deleted.is_(False)) if tree_ids: - stmt = stmt.where(TreeModel.id.in_(tree_ids)) + stmt = stmt.where(_id_in(TreeModel.id, tree_ids)) rows = (await self.db.execute(stmt)).all() if not rows: raise CustomException(msg="无活动单株") @@ -3203,11 +3220,35 @@ class StatisticsService: raise CustomException(msg=f"性状不可用于统计(需数值型且存在): {trait_code}") direction = (trait_row.direction or "desc").lower() desirable_high = direction != "asc" # desc/None:表型越高越好;asc:越低越好(如病指) - # 1) 单株系谱 → 深层系谱 A 矩阵 → 近交系数 F + # 0) 先取目标性状每树表型(同株多次取均值)——近交只针对有该性状表型的单株; + # 全库 3.6 万树建稠密 A 会 OOM(36990²×8B≈10 GiB),必须按表型树收口。 + obs = TraitObservationModel + stmt = ( + select(obs.tree_id, func.avg(obs.value_numeric).label("v")) + .join(TreeModel, obs.tree_id == TreeModel.id) + .where( + obs.trait_id == trait_id, + obs.value_numeric.isnot(None), + obs.tree_id.isnot(None), + obs.is_deleted.is_(False), + TreeModel.is_deleted.is_(False), + ) + .group_by(obs.tree_id) + ) + if year is not None: + stmt = stmt.where(obs.evaluate_year == year) + if trial_study_id is not None: + stmt = stmt.where(obs.trial_study_id == trial_study_id) + agg = (await self.db.execute(stmt)).all() + if not agg: + raise CustomException(msg="无可用表型数据(该性状为空)") + # 1) 单株系谱(仅表型树)→ 深层系谱 A 矩阵 → 近交系数 F # 单株 dam/sire 引用种质;种质系谱经 _germplasm_pedigree 多源 BFS 解析(bre_pedigree # 名称反查 + tree FK 兜底 + 祖先展开,与 OCS/mating 同源)——多世代项目(半同胞共享 # 祖亲)不再低估 F;被引用但未解析(软删等)种质补 founder 兜底。 - stmt = select(TreeModel.id, TreeModel.dam_id, TreeModel.sire_id).where(TreeModel.is_deleted.is_(False)) + stmt = (select(TreeModel.id, TreeModel.dam_id, TreeModel.sire_id) + .where(TreeModel.is_deleted.is_(False), + _id_in(TreeModel.id, [r.tree_id for r in agg]))) rows = (await self.db.execute(stmt)).all() if not rows: raise CustomException(msg="无活动单株") @@ -3249,27 +3290,7 @@ class StatisticsService: for k, v in rmat_res["inbreeding"].items(): if k.startswith("t") and k[1:].isdigit(): f_by_tree[int(k[1:])] = float(v) - # 2) 目标性状每树表型(同株多次取均值) - obs = TraitObservationModel - stmt = ( - select(obs.tree_id, func.avg(obs.value_numeric).label("v")) - .join(TreeModel, obs.tree_id == TreeModel.id) - .where( - obs.trait_id == trait_id, - obs.value_numeric.isnot(None), - obs.tree_id.isnot(None), - obs.is_deleted.is_(False), - TreeModel.is_deleted.is_(False), - ) - .group_by(obs.tree_id) - ) - if year is not None: - stmt = stmt.where(obs.evaluate_year == year) - if trial_study_id is not None: - stmt = stmt.where(obs.trial_study_id == trial_study_id) - agg = (await self.db.execute(stmt)).all() - if not agg: - raise CustomException(msg="无可用表型数据(该性状为空)") + # 2) 表型已在上方取好(agg),直接配对 F pairs: list[tuple[int, float, float]] = [] for r in agg: f = f_by_tree.get(int(r.tree_id)) @@ -3329,7 +3350,7 @@ class StatisticsService: interpretation = f"未检出近交衰退(F 对表型无{trend}影响)" tree_ids = [p[0] for p in pairs] trows = (await self.db.execute( - select(TreeModel.id, TreeModel.tree_no).where(TreeModel.id.in_(tree_ids)) + select(TreeModel.id, TreeModel.tree_no).where(_id_in(TreeModel.id, tree_ids)) )).all() tree_map = {tid: no for tid, no in trows} base.update({ @@ -3470,7 +3491,7 @@ class StatisticsService: if flagged_ids: for tid, tno in (await self.db.execute( select(TreeModel.id, TreeModel.tree_no) - .where(TreeModel.id.in_(flagged_ids)))).all(): + .where(_id_in(TreeModel.id, flagged_ids)))).all(): tree_no_map[int(tid)] = tno or "" excluded_ids = await self._excluded_tree_ids() mendelian_section = { @@ -3665,7 +3686,7 @@ class StatisticsService: gs = (await self.db.execute( select(BreedingGermplasmModel.id, BreedingGermplasmModel.cultivar_name, BreedingGermplasmModel.accession_no) - .where(BreedingGermplasmModel.id.in_(gids)) + .where(_id_in(BreedingGermplasmModel.id, gids)) )).all() cand_ids = [int(g.id) for g in gs] if not cand_ids: @@ -3793,7 +3814,7 @@ class StatisticsService: pvs = (await self.db.execute( select(PredictionValueModel).where( PredictionValueModel.prediction_id == ebv_pred.id, - PredictionValueModel.germplasm_id.in_(cand_ids), + _id_in(PredictionValueModel.germplasm_id, cand_ids), PredictionValueModel.predicted_value.isnot(None), ) )).scalars().all() @@ -3819,7 +3840,7 @@ class StatisticsService: mset = set(missing) off_tree = (await self.db.execute( select(TreeModel.id, TreeModel.dam_id, TreeModel.sire_id, TreeModel.combination_id).where( - or_(TreeModel.dam_id.in_(missing), TreeModel.sire_id.in_(missing)), + or_(_id_in(TreeModel.dam_id, missing), _id_in(TreeModel.sire_id, missing)), TreeModel.is_deleted.is_(False), ) )).all() @@ -3828,7 +3849,7 @@ class StatisticsService: off_pvs = (await self.db.execute( select(PredictionValueModel.tree_id, PredictionValueModel.predicted_value).where( PredictionValueModel.prediction_id == ebv_pred.id, - PredictionValueModel.tree_id.in_(off_ids), + _id_in(PredictionValueModel.tree_id, off_ids), PredictionValueModel.predicted_value.isnot(None), ) )).all() @@ -3887,7 +3908,7 @@ class StatisticsService: cand_ids = list(dict.fromkeys(candidate_germplasm_ids)) gs = (await self.db.execute( select(BreedingGermplasmModel).where( - BreedingGermplasmModel.id.in_(cand_ids), + _id_in(BreedingGermplasmModel.id, cand_ids), BreedingGermplasmModel.is_deleted.is_(False), ) )).scalars().all() @@ -4048,7 +4069,7 @@ class StatisticsService: raise CustomException(msg=f"选择数量 n_select={n_select} 超过候选亲本数 {len(cand_ids)}") gs = (await self.db.execute( select(BreedingGermplasmModel).where( - BreedingGermplasmModel.id.in_(cand_ids), + _id_in(BreedingGermplasmModel.id, cand_ids), BreedingGermplasmModel.is_deleted.is_(False), ) )).scalars().all() @@ -4210,7 +4231,7 @@ class StatisticsService: cand_rows = (await self.db.execute( select(TreeModel.id, TreeModel.tree_no).where( - TreeModel.id.in_(list(dict.fromkeys(candidate_tree_ids))), + _id_in(TreeModel.id, list(dict.fromkeys(candidate_tree_ids))), TreeModel.is_deleted.is_(False), ) )).all() @@ -4241,7 +4262,7 @@ class StatisticsService: select(GenotypeCallModel.sample_id, GenotypeCallModel.marker_id, GenotypeCallModel.allele) .where( - GenotypeCallModel.marker_id.in_(bg_marker_ids), + _id_in(GenotypeCallModel.marker_id, bg_marker_ids), GenotypeCallModel.allele.isnot(None), GenotypeCallModel.is_deleted.is_(False), ) @@ -4249,7 +4270,7 @@ class StatisticsService: if calls: sample_ids = sorted({c.sample_id for c in calls}) samples = (await self.db.execute( - select(GenotypeSampleModel).where(GenotypeSampleModel.id.in_(sample_ids)) + select(GenotypeSampleModel).where(_id_in(GenotypeSampleModel.id, sample_ids)) )).scalars().all() sample_tree, _ = await self._sample_to_tree_map(samples) for c in calls: @@ -4412,7 +4433,7 @@ class StatisticsService: block_of[ent.id] = blk await self.db.execute( update(TrialStudyEntryModel) - .where(TrialStudyEntryModel.id.in_(list(block_of))) + .where(_id_in(TrialStudyEntryModel.id, list(block_of))) .values(block_no=case(block_of, value=TrialStudyEntryModel.id, else_=None)) ) await self.db.flush() @@ -4468,7 +4489,7 @@ class StatisticsService: if block_of: await self.db.execute( update(TrialStudyEntryModel) - .where(TrialStudyEntryModel.id.in_(list(block_of))) + .where(_id_in(TrialStudyEntryModel.id, list(block_of))) .values(block_no=case(block_of, value=TrialStudyEntryModel.id, else_=None)) ) await self.db.flush() @@ -4577,7 +4598,7 @@ class StatisticsService: if block_of: await self.db.execute( update(TrialStudyEntryModel) - .where(TrialStudyEntryModel.id.in_(list(block_of))) + .where(_id_in(TrialStudyEntryModel.id, list(block_of))) .values(block_no=case(block_of, value=TrialStudyEntryModel.id, else_=None)) ) await self.db.flush() @@ -4634,7 +4655,7 @@ class StatisticsService: if block: tree_ids = [int(r.tree_id) for r in agg] block_rows = (await self.db.execute( - select(TreeModel.id, TreeModel.block_no).where(TreeModel.id.in_(tree_ids)) + select(TreeModel.id, TreeModel.block_no).where(_id_in(TreeModel.id, tree_ids)) )).all() block_of = {int(tid): bn for tid, bn in block_rows} rows = [ @@ -4998,10 +5019,10 @@ class StatisticsService: all_tree_ids = sorted(set(by_tree) | set(ebv_map) | mas["tree_ids"]) if all_tree_ids: tree_no_map = dict((await self.db.execute( - select(TreeModel.id, TreeModel.tree_no).where(TreeModel.id.in_(all_tree_ids)) + select(TreeModel.id, TreeModel.tree_no).where(_id_in(TreeModel.id, all_tree_ids)) )).all()) combo_map = dict((await self.db.execute( - select(TreeModel.id, TreeModel.combination_id).where(TreeModel.id.in_(all_tree_ids)) + select(TreeModel.id, TreeModel.combination_id).where(_id_in(TreeModel.id, all_tree_ids)) )).all()) else: tree_no_map, combo_map = {}, {} @@ -5168,7 +5189,7 @@ class StatisticsService: tree_ids = sorted({r.tree_id for r in rows}) trows = (await self.db.execute( select(TreeModel.id, TreeModel.clone_id, TreeModel.combination_id).where( - TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + _id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() clone_of = {t.id: t.clone_id for t in trows} combo_of = {t.id: t.combination_id for t in trows} @@ -5545,7 +5566,7 @@ class StatisticsService: else: tree_rows = (await self.db.execute( select(TreeModel.id, TreeModel.clone_id, TreeModel.tree_no).where( - TreeModel.id.in_(list(by_tree)), TreeModel.is_deleted.is_(False) + _id_in(TreeModel.id, list(by_tree)), TreeModel.is_deleted.is_(False) ) )).all() clone_of = {t.id: t.clone_id for t in tree_rows} @@ -5555,7 +5576,7 @@ class StatisticsService: if clone_ids: clones = (await self.db.execute( select(CloneModel.id, CloneModel.clone_code).where( - CloneModel.id.in_(clone_ids), CloneModel.is_deleted.is_(False) + _id_in(CloneModel.id, clone_ids), CloneModel.is_deleted.is_(False) ) )).all() clone_label = {c.id: c.clone_code for c in clones} @@ -5807,7 +5828,7 @@ class StatisticsService: rel = result["reliability"] ebv = result["ebv"] germ_of = dict((await self.db.execute( - select(TreeModel.id, TreeModel.germplasm_id).where(TreeModel.id.in_(tree_ids)) + select(TreeModel.id, TreeModel.germplasm_id).where(_id_in(TreeModel.id, tree_ids)) )).all()) direction = (trait_row.direction or "desc").lower() tree_ebv = [(tid, float(ebv.get(f"t{tid}", 0.0))) for tid in tree_ids] @@ -5898,7 +5919,7 @@ class StatisticsService: select(GenotypeCallModel.sample_id, GenotypeCallModel.marker_id, GenotypeCallModel.allele) .where( - GenotypeCallModel.sample_id.in_(sample_ids), + _id_in(GenotypeCallModel.sample_id, sample_ids), GenotypeCallModel.marker_id.isnot(None), GenotypeCallModel.allele.isnot(None), GenotypeCallModel.is_deleted.is_(False), @@ -5909,7 +5930,7 @@ class StatisticsService: marker_ids = sorted({r.marker_id for r in call_rows}) markers = (await self.db.execute( select(MarkerModel.id, MarkerModel.marker_name, MarkerModel.marker_type).where( - MarkerModel.id.in_(marker_ids), MarkerModel.is_deleted.is_(False)) + _id_in(MarkerModel.id, marker_ids), MarkerModel.is_deleted.is_(False)) )).all() if not markers: raise CustomException(msg="标记不存在") @@ -6068,7 +6089,7 @@ class StatisticsService: select(GenotypeCallModel.sample_id, GenotypeCallModel.marker_id, GenotypeCallModel.allele) .where( - GenotypeCallModel.sample_id.in_(sample_ids), + _id_in(GenotypeCallModel.sample_id, sample_ids), GenotypeCallModel.marker_id.isnot(None), GenotypeCallModel.allele.isnot(None), GenotypeCallModel.is_deleted.is_(False), @@ -6081,7 +6102,7 @@ class StatisticsService: select(MarkerModel.id, MarkerModel.marker_name, MarkerModel.chromosome, MarkerModel.position, MarkerModel.marker_type) - .where(MarkerModel.id.in_(marker_ids), MarkerModel.is_deleted.is_(False)) + .where(_id_in(MarkerModel.id, marker_ids), MarkerModel.is_deleted.is_(False)) )).all() if not markers: raise CustomException(msg="标记不存在") @@ -6594,7 +6615,7 @@ class StatisticsService: select(GenotypeCallModel.sample_id, GenotypeCallModel.marker_id, GenotypeCallModel.allele) .where( - GenotypeCallModel.marker_id.in_(marker_ids), + _id_in(GenotypeCallModel.marker_id, marker_ids), GenotypeCallModel.allele.isnot(None), GenotypeCallModel.is_deleted.is_(False), ) @@ -6603,7 +6624,7 @@ class StatisticsService: return out sample_ids = sorted({c.sample_id for c in calls}) samples = (await self.db.execute( - select(GenotypeSampleModel).where(GenotypeSampleModel.id.in_(sample_ids)) + select(GenotypeSampleModel).where(_id_in(GenotypeSampleModel.id, sample_ids)) )).scalars().all() sample_tree, _ = await self._sample_to_tree_map(samples) tree_geno: dict[int, dict[int, str]] = {} @@ -6744,7 +6765,7 @@ class StatisticsService: # GXE 分支:环境维度记录 → rec_phenos(与 run_ablup 同构) grp_rows = (await self.db.execute( select(TreeModel.id, TreeModel.clone_id, TreeModel.combination_id).where( - TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + _id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() clone_of = {t.id: t.clone_id for t in grp_rows} combo_of = {t.id: t.combination_id for t in grp_rows} @@ -6764,7 +6785,7 @@ class StatisticsService: else: site_rows = (await self.db.execute( select(TreeModel.id, TreeModel.trial_study_id).where( - TreeModel.id.in_(tree_ids), TreeModel.is_deleted.is_(False)) + _id_in(TreeModel.id, tree_ids), TreeModel.is_deleted.is_(False)) )).all() site_of = {t.id: t.trial_study_id for t in site_rows} rec_phenos = {} diff --git a/backend/app/api/v1/module_bre/tree/service.py b/backend/app/api/v1/module_bre/tree/service.py index 6f46327..13a7a0d 100644 --- a/backend/app/api/v1/module_bre/tree/service.py +++ b/backend/app/api/v1/module_bre/tree/service.py @@ -12,6 +12,7 @@ from app.utils.dict_util import DictLabelResolver, dict_value_to_label from app.utils.number_gen import NumberGenService from .crud import BreedingTreeCRUD +from .model import TreeModel from .schema import ( TreeCreateSchema, TreeOutSchema, @@ -370,8 +371,10 @@ class TreeService: await BreedingTreeCRUD(self.auth, self.db).delete(ids=ids) async def list_options(self) -> list[dict[str, Any]]: - """供前端下拉选择使用:返回 [{value, label}]。""" - obj_list = await BreedingTreeCRUD(self.auth, self.db).get_list(order_by=[{"id": "asc"}]) + """供前端下拉选择使用:返回 [{value, label}]。只取 id/tree_no 两列,避免全表列加载拖慢下拉。""" + obj_list = await BreedingTreeCRUD(self.auth, self.db).get_list( + order_by=[{"id": "asc"}], load_columns=[TreeModel.id, TreeModel.tree_no] + ) return [{"value": o.id, "label": o.tree_no} for o in obj_list] @staticmethod diff --git a/backend/scripts/breeding_stats/blup.py b/backend/scripts/breeding_stats/blup.py index 226525e..988cd9f 100644 --- a/backend/scripts/breeding_stats/blup.py +++ b/backend/scripts/breeding_stats/blup.py @@ -9,15 +9,21 @@ - 无父母的个体一律视为 founder(base population,彼此不相关)。 算法: -1. A 矩阵(稠密,O(n²) 递归)→ 各体自交系数 F; -2. A⁻¹ 用 Henderson 稀疏规则构造(无需求 A 的逆); -3. MME 求解:y = Xb + Zu + e,X=截距+固定效应虚拟列+协变量(默认仅截距),λ = σ²e/σ²a; - 默认稠密精确求解;个体数 > N_SPARSE(默认 1000)且系谱存在时切稀疏 A⁻¹ + 共轭梯度 - 迭代(solve 的 solver="sparse-cg"):EBV/h²/σ² 为 CG 精确解,reliability 用 - Hutchinson 随机探测估计 C22 对角线(reliability_approx=True,均值/PA 可靠、个体级近似); -4. 剖面 REML:对每个候选 h² 固定 λ,用 MME 解的 EM 条件式估计 σ²e、σ²a, - 再用精确 REML 对数似然评分,黄金分割一维最大化(稳健、可收敛到边界); -5. h² = σ²a/(σ²a+σ²e);reliability = 1 - PEV/σ²a,PEV = C22[ii]·σ²e。 +1. A⁻¹ 用 Henderson 稀疏规则构造(无需求 A 的逆); +2. MME 求解:y = Xb + Zu + e,X=截距+固定效应虚拟列+协变量(默认仅截距),λ = σ²e/σ²a; + 个体数 ≤ N_SPARSE(默认 1000)→ 稠密精确求解(solver="dense":索引化 Z、稠密 A、 + 剖面 REML 黄金分割);个体数 > N_SPARSE 且系谱存在 → 稀疏路径(solver="sparse-em"): + Z 索引化设计矩阵(bincount)、稀疏 A⁻¹ 三元组 + 共轭梯度迭代,全程不落稠密 Z/A/V + (35k 树在常规内存内跑通)。方差分量按 n 分档: + 个体数 ≤ N_EXACT(3000)→ tr(A⁻¹C²²) 逐列 CG 精确(无随机噪声),REML 解 + h(λ)=g(λ)−λ=0 单调单根 → log10λ 二分 ~32 步到 1e-5(规避 EM 慢收缩); + 个体数 > N_EXACT 且观测数 m>N_SUBSAMPLE(可真实子采样)→ 随机子样本(观测 + ≤N_SUBSAMPLE,含祖先闭包)上精确迹二分求 REML 解——子样本 REML 为全群方差组分 + 一致估计(抽样 SE h² 约 ±0.04),再全群 MME 于子样本 λ* 下解 EBV(大样本 EM 慢收缩 + ρ≈0.999 冻结不可用;EBV 对 λ 稳健,warning 标注); +3. reliability:稠密路径精确(diag(C22) 取 MME 逆对角);稀疏路径 n≤EXACT_DIAG_N 用 + 稠密 MME 逆精确、n 大用 Hutchinson 随机探测(reliability_approx=True,均值/PA 可靠、个体级近似); +4. h² = σ²a/(σ²a+σ²e);reliability = 1 - PEV/σ²a,PEV = C22[ii]·σ²e。 无系谱(所有个体均无父母)时:不背书遗传力—— h2 置 None、reliability 全部置 0、EBV 取中心化表型,并返回 warning。 @@ -41,8 +47,27 @@ CG_TOL = 1e-9 # 共轭梯度相对残差收敛阈值 # CG 精确解;reliability 是 MC 对角估计——同家族强相关的 C22 令其收敛慢(~1/√k), # 误差随 k 减小但非精确,作为文档标注的近似(排序可靠、量级可信)。 HUTCH_K = 100 +# 稀疏路径 reliability 分档:个体数 ≤ EXACT_DIAG_N → diag(C22) 用稠密 MME 逆精确求 +# (golden 等价档,无 MC 噪声——n=20 全同胞小样本 h²→1 边界实测 Hutchinson 均值误差 ~0.08); +# > EXACT_DIAG_N → Hutchinson 随机探测(排序可靠、量级可信,误差随个体数稀释)。 +EXACT_DIAG_N = 1000 +# 稀疏路径分档:精确迹二分求 λ*(逐列 CG 精确 tr(A⁻¹C²²),无随机噪声, +# REML 解 = h(λ)=g(λ)−λ=0 单调单根,log10λ 二分 ~32 步到 1e-5,规避 EM 慢收缩 ρ≈0.987)。 +# 精确迹 O(n²)(每步似然 n 次逐列 CG),仅个体数 ≤ N_EXACT 才可负担;>N_EXACT 且可真实 +# 子采样观测(m>N_SUBSAMPLE)→ 随机子样本精确迹二分(观测 ≤N_SUBSAMPLE,含祖先闭包) +# 求 λ*,再全群 MME 求解。 +# 大样本 EM 收缩率 ρ→0.999(100 步 λ 几乎不动,冻结在初始猜测),Hutchinson 噪声又淹没 +# h(λ) 信号(~2e-4 vs 噪声 ~2e-3),全群精确迹 O(n²) 不可行;子样本 iid 抽样 REML 为全群 +# 方差组分一致估计(n_sub=3000 抽样 SE h² ≈ ±0.04),EBV 对 λ 稳健(偏 25% Kendall τ>0.97)。 +# N_EXACT 对齐 N_SUBSAMPLE:m>N_SUBSAMPLE 时 n≈m+系谱必超 N_EXACT → 走子样本封顶 ~70min, +# 消除 3000N_EXACT 且 m>N_SUBSAMPLE 才触发) +BISECT_LO, BISECT_HI = -4.0, 4.0 # log10(λ=Ve/Va) 二分区间(REML 解通常在其内) +BISECT_TOL = 1e-8 # 二分 log10(λ) 收敛宽度(→ λ 相对精度 ~1e-8) +MAX_BISECT_EVAL = 60 # 二分求值上限 -ENGINE_VERSION = "1.5.0" # 单性状动物模型 BLUP 求解器版本(落 bre_prediction.engine_version;1.1.0: 供 mtblup 单性状方差复用;1.2.0: G×R 砧木×接穗互作;1.3.0: AR1×AR1 空间协方差 solve_spatial + 稀疏 A⁻¹/共轭梯度/Hutchinson 可靠性;1.4.0: AR1×AR1 各向异性双参数 ρ_row/ρ_col;1.5.0: 区组随机效应 _solve_block(增广/α-格子 block_no 第二随机效应)) +ENGINE_VERSION = "1.7.0" # 单性状动物模型 BLUP 求解器版本(落 bre_prediction.engine_version;1.1.0: 供 mtblup 单性状方差复用;1.2.0: G×R 砧木×接穗互作;1.3.0: AR1×AR1 空间协方差 solve_spatial + 稀疏 A⁻¹/共轭梯度/Hutchinson 可靠性;1.4.0: AR1×AR1 各向异性双参数 ρ_row/ρ_col;1.5.0: 区组随机效应 _solve_block(增广/α-格子 block_no 第二随机效应);1.6.0: 稀疏路径 EM-REML(索引化设计矩阵 + 稀疏 A⁻¹ 无稠密 A + EM-REML,35k 树不再 OOM;n≤N_EXACT 精确迹二分 λ*,n>N_EXACT 随机子样本精确迹二分 → 全群 MME 解 EBV);1.7.0: 性能倒挂修复——N_EXACT 6000→3000 对齐 N_SUBSAMPLE,生产档仅 m>N_SUBSAMPLE 才触发;系谱闭包改 BFS 可达祖先(服务层),选中 200~3000 树不再被全量 germplasm 闭包顶进数小时精确迹尾) def _pearson(x: list[float], y: list[float]) -> float: @@ -129,30 +154,37 @@ def _build_ainv(order: list[int], n_base: int, def _build_ainv_sparse(order: list[int], n_base: int, parent_of: dict[int, tuple[int | None, int | None]] ) -> tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]: - """返回 (A, ii, jj, vv):A 稠密(供 REML 似然与近交系数);A⁻¹ 以 COO 稀疏三元组存储。 + """返回 (diag_a, ii, jj, vv):A 对角(近交 F=A_ii−1 来源);A⁻¹ 以 COO 稀疏三元组存储。 - Henderson 规则与 _build_ainv 逐项一致,但 A⁻¹ 不落稠密矩阵,只收集非零元 - (对角 + 2 亲本 × (对角/双向共祖先)),供共轭梯度 matvec 使用: - (Ainv·u)[i] = Σ_j Ainv_ij·u_j,np.add.at 按三元组累加即可。 + 不构建稠密 A:加性关系对角用「祖先稀疏行」递推(每体一行 {祖先: A_ij},浅系谱祖先数 + 小,O(n·s²)),Henderson 规则只需 A 对角(fp/fq)与双亲交叉 A[di,si],均可在稀疏行上 + 读出;A⁻¹ 非零元(对角 + 2 亲本 × 对角/双向共祖先)逐项与 _build_ainv 一致。 + matvec 时用 np.bincount 按 ii 累加 ainv_vv·u[jj] 即得 A⁻¹·u(重复元自动合并)。 """ n = len(order) idx = {ind: pos for pos, ind in enumerate(order)} - A = np.eye(n) + anc: list[dict[int, float]] = [None] * n + diag_a = np.zeros(n) for pos, ind in enumerate(order): d, s = parent_of.get(ind, (None, None)) if d is None and s is None: + anc[pos] = {pos: 1.0} + diag_a[pos] = 1.0 continue di = idx.get(d) if d is not None else None si = idx.get(s) if s is not None else None - row = np.zeros(n) - if di is not None: - row += A[di] - if si is not None: - row += A[si] - A[pos, :] = 0.5 * row - A[:pos, pos] = A[pos, :pos] - apq = A[di, si] if (di is not None and si is not None) else 0.0 - A[pos, pos] = 1.0 + 0.5 * apq + row: dict[int, float] = {} + for src in (di, si): + if src is not None: + for j, v in anc[src].items(): + row[j] = row.get(j, 0.0) + 0.5 * v + apq = 0.0 + if di is not None and si is not None: + small, big = (anc[di], anc[si]) if len(anc[di]) < len(anc[si]) else (anc[si], anc[di]) + apq = sum(v * big[j] for j, v in small.items() if j in big) + row[pos] = 1.0 + 0.5 * apq + anc[pos] = row + diag_a[pos] = row[pos] ii: list[int] = [] jj: list[int] = [] @@ -171,8 +203,8 @@ def _build_ainv_sparse(order: list[int], n_base: int, i = idx[ind] di = idx.get(d) if d is not None else None si = idx.get(s) if s is not None else None - fp = A[di, di] - 1.0 if di is not None else 0.0 - fq = A[si, si] - 1.0 if si is not None else 0.0 + fp = diag_a[di] - 1.0 if di is not None else 0.0 + fq = diag_a[si] - 1.0 if si is not None else 0.0 t = 0.5 - 0.25 * (fp + fq) # 极端近交(双亲 F→1)时 t→0,除零保护;1e-6 下界对正常系谱无感 inv_t = 1.0 / max(t, 1e-6) @@ -185,7 +217,7 @@ def _build_ainv_sparse(order: list[int], n_base: int, if di is not None and si is not None: _put(di, si, 0.25 * inv_t) _put(si, di, 0.25 * inv_t) - return (A, np.array(ii, dtype=np.int64), np.array(jj, dtype=np.int64), + return (diag_a, np.array(ii, dtype=np.int64), np.array(jj, dtype=np.int64), np.array(vv, dtype=float)) @@ -351,11 +383,13 @@ def solve(pedigree: list[dict], phenotypes: dict, *, tol: float = TOL, gxe_ratio / n_gxe / n_cross_env / n_genotypes(G×E 分支); block_ratio / n_blocks / block_effects(区组分支); n_obs / n_individuals / n_base / n_with_parents / n_fixed: int; - converged: bool;n_iter: int(剖面求值次数/EM 迭代数);warning: str | None; - solver: str("dense" 稠密精确 / "sparse-cg" 稀疏共轭梯度,n>N_SPARSE 时触发); - 稀疏路径另含 cg_iter(CG 迭代数)与 reliability_approx=True(可靠性为 Hutchinson 近似, - EBV/h²/σ² 仍为精确解)。 + converged: bool;n_iter: int(稠密=剖面求值次数,稀疏=迹二分/子样本求值次数);warning: str | None; + solver: str("dense" 稠密精确 / "sparse-em" 稀疏 EM-REML,n>N_SPARSE 时触发); + 稀疏路径另含 cg_iter(CG 迭代数)与 reliability_approx=True(可靠性为 Hutchinson 近似); + 生产档(n>N_EXACT)方差组分为随机子样本 REML 一致估计(抽样 SE h² ≈ ±0.04),EBV 由 + 全群 MME 于子样本 λ 下求解——对 λ 稳健,排序可靠,h²/σ² 视为近似区间。 """ + global N_EXACT # 生产档内联递归临时抬档(try/finally 还原);须在函数内首次读取 N_EXACT 前声明 if gxe is not None: return _solve_gxe(pedigree, phenotypes, fixed=fixed, covariate=covariate, gxe=gxe, record_map=record_map or {}, tol=tol, @@ -427,22 +461,237 @@ def solve(pedigree: list[dict], phenotypes: dict, *, tol: float = TOL, "warning": "无系谱信息(所有个体均无父母):EBV 仅为表型残差(固定效应已校正),未估计遗传力与可靠性。", } + XtX = X.T @ X + Xty = X.T @ y + # 索引化设计矩阵(Z 每行恰一个 1):Zty / ZtZ 对角 / XtZ 全用 bincount 累加, + # 稀疏路径不再落稠密 Z(35k 树稠密 Z ≈ 54 GiB)。稠密分支在小 n 下重建等价 Z。 + obs_idx = np.array([idx[ind] for ind in obs_ids], dtype=np.int64) + Zty = np.bincount(obs_idx, weights=y, minlength=n) + ZtZ_diag = np.bincount(obs_idx, minlength=n).astype(float) + XtZ = np.zeros((p, n)) + for col in range(p): + XtZ[col] = np.bincount(obs_idx, weights=X[:, col], minlength=n) + + # A7 阈值分派:个体数 > N_SPARSE 且系谱存在 → 稀疏 A⁻¹ + 索引化设计 + EM-REML + # (回避 O(n³) 稠密 MME 逆与 O(m²) 稠密 V);可靠性 n≤EXACT_DIAG_N 用稠密 MME 逆精确, + # n 大用 Hutchinson 随机探测估计 C22 对角线。 + use_sparse = n_with_parents > 0 and n > N_SPARSE + if use_sparse: + diag_a, ainv_ii, ainv_jj, ainv_vv = _build_ainv_sparse(order, n_base, parent_of) + else: + A, Ainv = _build_ainv(order, n_base, parent_of) + + if use_sparse: + # ---- 稀疏路径:EM-REML(Meyer 1985 单性状动物模型)---- + # Va = (û'A⁻¹û + σ²e·tr(A⁻¹C²²))/n;Ve = yPy/(m−p);λ=Ve/Va。全程不落稠密 Z/A/V。 + # n ≤ N_EXACT:tr(A⁻¹C²²) 逐列 CG 精确(无随机噪声);REML 解 = h(λ)=g(λ)−λ=0 + # 单调单根 → log10λ 二分,~32 步到 1e-5(规避 EM 慢收缩 ρ≈0.987); + # n > N_EXACT:随机子样本(≤N_SUBSAMPLE 观测,祖先闭包)内联递归走精确迹二分 + # 求 λ*,再全群 MME 于 λ* 解 EBV——子样本抽样误差即真实不确定性(warning 标注)。 + ztz = ZtZ_diag + + def _matvec(lam: float): + def mv(v: np.ndarray) -> np.ndarray: + b, u = v[:p], v[p:] + out1 = XtX @ b + XtZ @ u + out2 = XtZ.T @ b + ztz * u + # A⁻¹·u:bincount 按 ii 累加(重复元自动合并),比 np.add.at 快一个量级 + out2 += np.bincount(ainv_ii, weights=lam * ainv_vv * u[ainv_jj], minlength=n) + return np.concatenate([out1, out2]) + return mv + + mme_rhs = np.concatenate([Xty, Zty]) + rng = np.random.RandomState(20260804) + + em_it = 0 + em_ok = True + cg_ok = True + cg_iter = 0 + # 生产档仅当能真实子采样观测(m>N_SUBSAMPLE)才走:m 不大时子样本=全样本,退化为 + # 精确迹 + 冗余全群 MME(更慢且子样本误差 SE 无意义)——直接全 n 精确迹更省更准。 + use_production = n > N_EXACT and m > N_SUBSAMPLE + if not use_production: + # ---- 精确迹二分(golden 等价档)---- + # A⁻¹ 第 j 列 = {i: ainv[i,j]};C²²e_j = MME 解 RHS=[0;e_j] 的 u 块; + # tr(A⁻¹C²²) = Σ_j (A⁻¹e_j)ᵀ(C²²e_j)。逐列 CG 精确,无 Hutchinson 噪声。 + col_of: list[list[tuple[int, float]]] = [[] for _ in range(n)] + for a, b, v in zip(ainv_ii.tolist(), ainv_jj.tolist(), ainv_vv.tolist()): + col_of[b].append((a, v)) + + def _em_step(lam: float) -> tuple[float, float]: + mv = _matvec(lam) + sol, _ci, _ok = _cg_solve(mv, mme_rhs) + b, u = sol[:p], sol[p:] + yPy = float(y @ y - b @ Xty - u @ Zty) + Ve = max(yPy / max(m - p, 1), _FLOOR) + tr = 0.0 + for j in range(n): + rhs2 = np.zeros(p + n) + rhs2[p + j] = 1.0 + sol2, _i2, _o2 = _cg_solve(mv, rhs2) + uj = sol2[p:] + s = 0.0 + for i, v in col_of[j]: + s += v * uj[i] + tr += s + uAu = float(np.bincount(ainv_ii, weights=ainv_vv * u[ainv_jj], minlength=n) @ u) + Va = max((uAu + Ve * tr) / max(n, 1), _FLOOR) + return Va, Ve + + def _h(lv: float) -> float: + lam = 10.0 ** lv + Va, Ve = _em_step(lam) + return Ve / Va - lam + + lo, hi = BISECT_LO, BISECT_HI + h_lo, h_hi = _h(lo), _h(hi) + while h_lo <= 0.0 and lo > BISECT_LO - 4.0: + lo -= 2.0 + h_lo = _h(lo) + while h_hi >= 0.0 and hi < BISECT_HI + 4.0: + hi += 2.0 + h_hi = _h(hi) + if h_lo > 0.0 > h_hi: + em_it = 2 + while (hi - lo) > BISECT_TOL and em_it < MAX_BISECT_EVAL: + mid = 0.5 * (lo + hi) + h_mid = _h(mid) + em_it += 1 + if h_mid > 0.0: + lo = mid + else: + hi = mid + Va, Ve = _em_step(10.0 ** (0.5 * (lo + hi))) + else: + # 无变号(REML 解在区间外/边界最优):取 |h| 最小的网格点 + lv_grid = np.linspace(max(lo - 4.0, -8.0), min(hi + 4.0, 8.0), 41) + h_abs = [abs(_h(lv)) for lv in lv_grid] + Va, Ve = _em_step(10.0 ** float(lv_grid[int(np.argmin(h_abs))])) + em_ok = False + else: + # ---- 生产档(n>N_EXACT 且 m>N_SUBSAMPLE):随机子样本精确迹二分 → 全群 MME ---- + # 全群精确迹 O(n²) 不可行;EM 收缩率 ρ→0.999 冻结在初始猜测,Hutchinson 噪声 + # 又淹没 h(λ) 信号(~2e-4 vs 噪声 ~2e-3)。子样本 iid 抽样 REML 为全群方差组分 + # 一致估计(n_sub=3000 抽样 SE h² ≈ ±0.04,EBV corr>0.998),抽样误差即真实 + # 不确定性。内联递归把 N_EXACT 临时抬到 10⁹,让子样本走精确迹二分求 λ*。 + chosen = set(rng.choice(obs_ids, size=min(N_SUBSAMPLE, m), replace=False).tolist()) + # 祖先闭包:parent_of 以个体 id 为键、值亦为 id(与 _build_ainv_sparse 一致); + # 勿用 idx[]/order[] 做中间层(那是位置空间,混用会把闭包加进错对象)。 + q: list[int] = list(chosen) + while q: + cur = q.pop() + d, s = parent_of.get(cur, (None, None)) + for pp in (d, s): + if pp is not None and pp not in chosen: + chosen.add(pp) + q.append(pp) + ped_sub = [] + for i in chosen: + d, s = parent_of.get(i, (None, None)) + ped_sub.append({ + "individual": i, + "dam": d if d is not None else None, + "sire": s if s is not None else None, + }) + phen_sub = {k: v for k, v in phenotypes.items() if k in chosen} + fixed_sub = None + if fixed: + fixed_sub = {fk: {k: v for k, v in fv.items() if k in chosen} + for fk, fv in fixed.items()} + cov_sub = None + if covariate: + cov_sub = {k: v for k, v in covariate.items() if k in chosen} + saved_exact = N_EXACT + try: + N_EXACT = 10 ** 9 # noqa: PLW0603 临时抬档(子样本个体数远小于全群) + sub = solve(ped_sub, phen_sub, fixed=fixed_sub, covariate=cov_sub) + finally: + N_EXACT = saved_exact + Va = max(float(sub["sigma_a"]), _FLOOR) + Ve = max(float(sub["sigma_e"]), _FLOOR) + em_ok = bool(sub.get("converged", False)) + em_it = int(sub.get("n_iter", 0)) + + lam = Ve / Va + mv = _matvec(lam) + sol, cg_iter, cg_ok = _cg_solve(mv, mme_rhs) + u = sol[p:] + h2 = Va / (Va + Ve) if (Va + Ve) > 0 else None + + # 可靠性:n ≤ EXACT_DIAG_N(稠密逆仍便宜)→ diag(C22) 用稠密 MME 逆精确求(golden + # 等价档,无 MC 噪声——全同胞小样本 h²→1 边界 Hutchinson 均值误差实测 ~0.08 量级, + # 见 e2e_spatial_sparse n=20 段);n 大时退 Hutchinson 随机探测(排序可靠、量级可信)。 + if n <= EXACT_DIAG_N: + Ainv_dense = np.zeros((n, n)) + np.add.at(Ainv_dense, (ainv_ii, ainv_jj), ainv_vv) # 三元组含重复 (i,j)(多后代),须累加 + big = np.zeros((p + n, p + n)) + big[:p, :p] = XtX + big[:p, p:] = XtZ + big[p:, :p] = XtZ.T + big[p:, p:] = np.diag(ZtZ_diag) + lam * Ainv_dense + try: + big_inv = np.linalg.inv(big) + except np.linalg.LinAlgError: + big_inv = np.linalg.pinv(big) + C22d = big_inv[p:, p:] + pev = np.diag(C22d) * Ve + else: + diag_c22 = np.zeros(n) + for _ in range(HUTCH_K): + pv = rng.choice([-1.0, 1.0], size=n) + rhs2 = np.zeros(p + n) + rhs2[p:] = pv + sol2, _it2, _ok2 = _cg_solve(mv, rhs2) + diag_c22 += sol2[p:] * pv + diag_c22 /= HUTCH_K + pev = diag_c22 * Ve + + ebv = {ind: float(u[idx[ind]]) for ind in order} + rel = {ind: float(np.clip(1.0 - pev[idx[ind]] / Va, 0.0, 1.0)) for ind in order} + + warning = None + if use_production: + # 生产档(n>N_EXACT 且 m>N_SUBSAMPLE):方差组分来自随机子样本 + # (≤N_SUBSAMPLE 观测,含祖先闭包)的精确 REML 估计——一致估计但带抽样误差 + # (n_sub=3000 SE h² ≈ ±0.04);EBV 由全群 MME 在子样本 λ* 下求解, + # 对 λ 稳健(偏 25% 内 Kendall τ>0.97),排序可靠。 + warning = (f"生产档(个体数 {n})方差分量取自随机子样本" + f"(≤{N_SUBSAMPLE} 观测,含祖先闭包)的 REML 一致估计:存在抽样误差" + "(h² 标准误约 ±0.04),h²/σ² 报告值请视为近似区间;EBV 由全群 MME" + "在子样本 λ 下求解,排名可靠。") + elif not em_ok: + warning = ("REML 解落在 log10λ 二分区间外(似然面边界最优)," + "结果采用 |h| 最小网格点,可能不精确。") + if not cg_ok: + resid = mme_rhs - mv(sol) + rel_res = float(np.linalg.norm(resid) / (np.linalg.norm(mme_rhs) + 1e-12)) + w = (f"共轭梯度 {cg_iter} 次迭代未达收敛阈值(相对残差 {rel_res:.2e})," + "结果采用已探明最优,可靠性为 Hutchinson 近似。") + warning = f"{warning} {w}" if warning else w + return { + "ebv": ebv, + "reliability": rel, + "h2": float(h2) if h2 is not None else None, + "sigma_a": float(Va), + "sigma_e": float(Ve), + "n_obs": m, + "n_individuals": n, + "n_base": n_base, + "n_with_parents": n_with_parents, + "n_fixed": n_fixed, + "converged": em_ok and cg_ok, + "n_iter": em_it, + "solver": "sparse-em", + "cg_iter": cg_iter, + "reliability_approx": n > EXACT_DIAG_N, + "warning": warning, + } + + # 稠密分支(n ≤ N_SPARSE):重建等价 Z / ZtZ(bincount 已给各体观测计数) Z = np.zeros((m, n)) for k, ind in enumerate(obs_ids): Z[k, idx[ind]] = 1.0 - XtX = X.T @ X - Xty = X.T @ y - Zty = Z.T @ y - ZtZ = Z.T @ Z - XtZ = X.T @ Z - - # A7 阈值分派:个体数 > N_SPARSE 且系谱存在 → 稀疏 A⁻¹ + 共轭梯度迭代求解 - # (回避 O(n³) 稠密 MME 逆);可靠性用 Hutchinson 随机探测估计 C22 对角线。 - use_sparse = n_with_parents > 0 and n > N_SPARSE - if use_sparse: - A, ainv_ii, ainv_jj, ainv_vv = _build_ainv_sparse(order, n_base, parent_of) - else: - A, Ainv = _build_ainv(order, n_base, parent_of) + ZtZ = np.diag(ZtZ_diag) def _reml_ll(Va: float, Ve: float) -> float: """精确 REML 对数似然(直接 V 计算,O(m³))。""" @@ -456,85 +705,6 @@ def solve(pedigree: list[dict], phenotypes: dict, *, tol: float = TOL, _, lx = np.linalg.slogdet(XtVinvX) return -0.5 * (float(lv) + float(lx) + yPy) - if use_sparse: - def _matvec(lam: float): - def mv(v: np.ndarray) -> np.ndarray: - b, u = v[:p], v[p:] - out1 = XtX @ b + XtZ @ u - out2 = XtZ.T @ b + ZtZ @ u - np.add.at(out2, ainv_ii, lam * ainv_vv * u[ainv_jj]) - return np.concatenate([out1, out2]) - return mv - - mme_rhs = np.concatenate([Xty, Zty]) - - best = {"h2": None, "ll": -np.inf, "Va": 0.0, "Ve": 0.0} - - def _profile(h2: float) -> float: - lam = (1.0 - h2) / h2 - sol, _it, _ok = _cg_solve(_matvec(lam), mme_rhs) - b = sol[:p] - u = sol[p:] - yPy = float(y @ y - b @ Xty - u @ Zty) - Ve = max(yPy / max(m - p, 1), _FLOOR) - Va = Ve / lam - ll = _reml_ll(Va, Ve) - if ll > best["ll"]: - best.update(h2=h2, ll=ll, Va=Va, Ve=Ve) - return ll - - h2_opt, n_iter = _golden_max(_profile, H2_MIN, H2_MAX, tol=tol) - _profile(h2_opt) - h2 = best["h2"] if best["h2"] is not None else h2_opt - Va, Ve = best["Va"], best["Ve"] - - lam = (1.0 - h2) / h2 - mv = _matvec(lam) - sol, cg_iter, cg_ok = _cg_solve(mv, mme_rhs) - u = sol[p:] - # Hutchinson:E[P⊙(C22·P)] = diag(C22),C22·P 由 CG 解 big·x=[0;P] 取 x[p:] 得到 - rng = np.random.RandomState(20260804) - diag_c22 = np.zeros(n) - for _ in range(HUTCH_K): - pv = rng.choice([-1.0, 1.0], size=n) - rhs2 = np.zeros(p + n) - rhs2[p:] = pv - sol2, _it2, _ok2 = _cg_solve(mv, rhs2) - diag_c22 += sol2[p:] * pv - diag_c22 /= HUTCH_K - - ebv = {ind: float(u[idx[ind]]) for ind in order} - pev = diag_c22 * Ve - rel = {ind: float(np.clip(1.0 - pev[idx[ind]] / Va, 0.0, 1.0)) for ind in order} - - warning = None - if n_iter >= MAX_PROFILE_EVALS: - warning = "剖面 REML 未完全收敛(似然面极平/边界最优),结果已采用已探明最优。" - if not cg_ok: - resid = mme_rhs - mv(sol) - rel_res = float(np.linalg.norm(resid) / (np.linalg.norm(mme_rhs) + 1e-12)) - w = (f"共轭梯度 {cg_iter} 次迭代未达收敛阈值(相对残差 {rel_res:.2e})," - "结果采用已探明最优,可靠性为 Hutchinson 近似。") - warning = f"{warning} {w}" if warning else w - return { - "ebv": ebv, - "reliability": rel, - "h2": float(h2), - "sigma_a": float(Va), - "sigma_e": float(Ve), - "n_obs": m, - "n_individuals": n, - "n_base": n_base, - "n_with_parents": n_with_parents, - "n_fixed": n_fixed, - "converged": n_iter < MAX_PROFILE_EVALS and cg_ok, - "n_iter": n_iter, - "solver": "sparse-cg", - "cg_iter": cg_iter, - "reliability_approx": True, - "warning": warning, - } - def _mme_solve(lam: float) -> tuple[np.ndarray, np.ndarray]: big = np.zeros((p + n, p + n)) big[:p, :p] = XtX diff --git a/backend/scripts/gen_obs_doc.py b/backend/scripts/gen_obs_doc.py new file mode 100644 index 0000000..df48d47 --- /dev/null +++ b/backend/scripts/gen_obs_doc.py @@ -0,0 +1,275 @@ +# -*- coding: utf-8 -*- +"""生成《桃育种观测业务梳理》Word 文档(微软雅黑,含表格与页码)。 + +运行:C:/ai/miniconda3/envs/dpb/python.exe backend/scripts/gen_obs_doc.py +输出:doc/桃育种观测业务梳理.docx +""" +import os +from docx import Document +from docx.shared import Pt, Cm, RGBColor +from docx.enum.text import WD_ALIGN_PARAGRAPH +from docx.enum.table import WD_TABLE_ALIGNMENT +from docx.oxml import OxmlElement +from docx.oxml.ns import qn + +FONT = "微软雅黑" +DARK = RGBColor(0x1F, 0x3B, 0x5C) +ROOT = os.path.dirname(os.path.dirname(os.path.dirname(os.path.abspath(__file__)))) +OUT = os.path.join(ROOT, "doc", "桃育种观测业务梳理.docx") + + +def set_run(run, size=10.5, bold=False, color=None): + run.font.name = FONT + run.font.size = Pt(size) + run.font.bold = bold + if color is not None: + run.font.color.rgb = color + run._element.rPr.rFonts.set(qn("w:eastAsia"), FONT) + + +def setup_styles(doc): + normal = doc.styles["Normal"] + normal.font.name = FONT + normal.font.size = Pt(10.5) + normal._element.rPr.rFonts.set(qn("w:eastAsia"), FONT) + for i in range(1, 7): + st = doc.styles["Heading %d" % i] + st.font.name = FONT + st.font.bold = True + st.font.color.rgb = DARK + st.font.size = Pt(max(16 - i * 1.5, 11)) + st._element.rPr.rFonts.set(qn("w:eastAsia"), FONT) + + +def add_title(doc, text, subtitle): + p = doc.add_paragraph() + p.alignment = WD_ALIGN_PARAGRAPH.CENTER + r = p.add_run(text) + set_run(r, size=20, bold=True, color=DARK) + p.paragraph_format.space_after = Pt(4) + p2 = doc.add_paragraph() + p2.alignment = WD_ALIGN_PARAGRAPH.CENTER + r2 = p2.add_run(subtitle) + set_run(r2, size=11, color=RGBColor(0x60, 0x60, 0x60)) + p2.paragraph_format.space_after = Pt(14) + + +def heading(doc, text, level): + p = doc.add_heading(text, level=level) + for r in p.runs: + set_run(r, bold=True) + p.paragraph_format.space_before = Pt(10 if level <= 2 else 6) + p.paragraph_format.space_after = Pt(4) + return p + + +def para(doc, text, size=10.5, bold=False, style=None, space_after=6): + p = doc.add_paragraph(style=style) + if text: + r = p.add_run(text) + set_run(r, size=size, bold=bold) + p.paragraph_format.space_after = Pt(space_after) + p.paragraph_format.line_spacing = 1.35 + return p + + +def shade_cell(cell, hexcolor): + tcPr = cell._tc.get_or_add_tcPr() + shd = OxmlElement("w:shd") + shd.set(qn("w:val"), "clear") + shd.set(qn("w:color"), "auto") + shd.set(qn("w:fill"), hexcolor) + tcPr.append(shd) + + +def table(doc, headers, rows, col_widths=None, size=9.5, header_fill="DCE6F1"): + t = doc.add_table(rows=1, cols=len(headers)) + t.style = "Table Grid" + t.alignment = WD_TABLE_ALIGNMENT.CENTER + for i, h in enumerate(headers): + cell = t.rows[0].cells[i] + p = cell.paragraphs[0] + r = p.add_run(h) + set_run(r, size=size, bold=True) + p.paragraph_format.space_after = Pt(0) + shade_cell(cell, header_fill) + for row in rows: + cells = t.add_row().cells + for i, v in enumerate(row): + p = cells[i].paragraphs[0] + r = p.add_run(str(v)) + set_run(r, size=size) + p.paragraph_format.space_after = Pt(0) + if col_widths: + for i, w in enumerate(col_widths): + for row in t.rows: + row.cells[i].width = Cm(w) + sp = doc.add_paragraph() + sp.paragraph_format.space_after = Pt(2) + return t + + +def add_page_number(section): + footer = section.footer + p = footer.paragraphs[0] + p.alignment = WD_ALIGN_PARAGRAPH.CENTER + r = p.add_run() + f1 = OxmlElement("w:fldChar"); f1.set(qn("w:fldCharType"), "begin") + it = OxmlElement("w:instrText"); it.set(qn("xml:space"), "preserve"); it.text = "PAGE" + f2 = OxmlElement("w:fldChar"); f2.set(qn("w:fldCharType"), "end") + r._r.append(f1); r._r.append(it); r._r.append(f2) + set_run(r, size=9, color=RGBColor(0x80, 0x80, 0x80)) + + +def main(): + doc = Document() + setup_styles(doc) + sec = doc.sections[0] + sec.page_width, sec.page_height = Cm(21.0), Cm(29.7) + sec.left_margin = sec.right_margin = Cm(2.2) + sec.top_margin = sec.bottom_margin = Cm(2.2) + add_page_number(sec) + + add_title(doc, "桃育种系统 · 育种观测业务梳理", "2026-08-06 · 内部梳理稿 v1.0") + + # 1 背景与目的 + heading(doc, "1 背景与目的", 1) + para(doc, "本次梳理源于一个具体问题:单株评价的新增/编辑表单中,「基本信息」之外的六个标签页(基本鉴定、果实外观、果皮、果实大小、果肉、果核)全部显示“暂无性状”,无法输入和修改。") + para(doc, "排查后发现这并非单个页面 bug,而是性状字典、观测数据模型与前端消费三者之间系统性不一致的暴露:库里现存 20 条统计性状的分类(产量/品质/果形/营养/抗性/物候),与表单六个标签页所用的官方描述字典(43 条,六大类)不匹配,导致标签页空转。") + para(doc, "本文档把围绕“观测”的现状梳理清楚:数据模型、核心表关系、三类观测切片、两条繁殖路径、扩繁追溯、性能、性状字典现状与断裂、已知缺口与待定决策,为后续方案设计提供一份完整、可复用的依据。") + + # 2 长表设计 + heading(doc, "2 观测数据模型:长表设计", 1) + para(doc, "核心观测表为 bre_trait_observation,采用“长表”(tidy)设计:每(单株 × 性状 × 年份)一行,一个性状一条记录。") + table(doc, + ["对比维度", "长表(当前)", "宽表(一性状一列)"], + [["数据结构", "每行一个观测值 + trait_id 外键", "每个性状一个列"], + ["性状可扩展", "加性状 = 加字典行,无需改表", "加性状 = 加列,需 DDL 迁移"], + ["页面/统计解耦", "页面按 trait 过滤即可", "列结构绑定页面"], + ["硬编码风险", "低(性状走字典)", "高"], + ["记录规模", "单表 358 万", "同规模列更多"]], + col_widths=[3.5, 6.5, 6.5]) + para(doc, "为什么必须长表:育种性状集合逐年演进(新性状、新标准、分子性状),宽表每次演进都要迁移;长表加一行字典即可,统计引擎按 data_type=='numeric' 直接消费。代价是记录量大,但见第 6 节,性能可用索引兜住。") + + # 3 四张核心表 + heading(doc, "3 四张核心表的关系", 1) + heading(doc, "3.1 单株评价 bre_tree_evaluation(主表)", 2) + para(doc, "一棵单株在某一年的一次评价(“头”)。记录评价年份、总评、负载量、评价人等信息。一个评价头对应一批性状观测明细。") + heading(doc, "3.2 性状观测 bre_trait_observation(明细表)", 2) + para(doc, "评价的“细目”,全系统唯一承载 358 万条观测的表。每行一个(单株 × 性状 × 年份)观测值,三态值列并存:") + para(doc, "· value_numeric —— 数值型(消费面:统计引擎)", style="List Bullet") + para(doc, "· value_text —— 文本 / 等级描述", style="List Bullet") + para(doc, "· value_date —— 日期型(物候等)", style="List Bullet") + heading(doc, "3.3 通用观测 bre_observation(独立轻量记录)", 2) + para(doc, "与“单株评价”解耦的散点式观测(园区环境、非评价场景的随记),trait 可空,不进主细结构。") + heading(doc, "3.4 性状字典 bre_trait", 2) + para(doc, "统一性状字典:code、中文名、unit、category、data_type、stage、方向(升/降)、是否入 EBV、h²、权重等。消费面是 data_type=='numeric' 的统计引擎。") + heading(doc, "3.5 主细关系", 2) + table(doc, + ["要点", "约定"], + [["基数", "bre_tree_evaluation 1 —— N bre_trait_observation"], + ["外键", "evaluation_id,CASCADE(删评价头连带删明细)"], + ["唯一约束", "UNIQUE(evaluation_id, trait_id):一个评价下同一性状只允许一条"]], + col_widths=[3.5, 13.0]) + + # 4 三类切片 + heading(doc, "4 观测数据的三类业务切片", 1) + para(doc, "观测行按「是否有评价头(evaluation_id)× stage」切为三类,合计 3,583,584 条。") + table(doc, + ["切片", "判定", "语义", "记录数", "备注"], + [["A 随评", "evaluation_id 非空", "随单株评价录入的性状观测", "731,904", "45,744 个评价头 × 16 条/头,0 孤儿"], + ["B 童期", "evaluation_id IS NULL + stage='juvenile'", "童期幼树的观测", "1,148,160", "无评价头,独立批次"], + ["C 克隆苗", "evaluation_id IS NULL + stage='evaluation'", "扩繁产出无性系苗的观测", "1,703,520", "每克隆 4 棵 ramet"]], + col_widths=[1.8, 5.0, 4.4, 2.2, 4.0]) + para(doc, "A 与 B/C 的判定只差一个字段是否为空,是历史演进留下的两类形态,方案阶段需决策是否统一(见第 9 节决策 3)。") + + # 5 育种流程全景 + heading(doc, "5 育种流程全景", 1) + heading(doc, "5.1 有性繁殖路径(实生苗)", 2) + para(doc, "杂交(cross_combination)→ 采粉(pollen)→ 授粉(pollination)→ 种子批(seed_lot)→ 种子处理(seed_treatment)→ 育苗(seedling)→ 种植(planting)→ 单株(tree)。种子苗是实生后代,遗传上每株独一无二,是选育的主体。") + heading(doc, "5.2 无性繁殖路径(克隆苗)", 2) + para(doc, "入选单株(selection_result)→ 晋级为克隆(tree_generation='clone')→ 扩繁(propagation,嫁接/扦插)→ 产出无性系苗 ramet(每克隆 4 棵)→ 克隆苗观测。克隆苗与母株基因型相同,用于性状复测、对比试验、示范与推广。") + heading(doc, "5.3 育苗 vs 克隆", 2) + table(doc, + ["维度", "育苗(有性)", "克隆扩繁(无性)"], + [["后代来源", "种子批 seed_lot", "接穗 scion_source(原株)"], + ["遗传构成", "实生后代,各不相同", "与母株基因型相同"], + ["用途", "选育新材料", "复测 / 示范 / 推广"], + ["观测归属", "单株 / 童期", "C 类克隆苗观测"]], + col_widths=[3.5, 6.5, 6.5]) + heading(doc, "5.4 扩繁追溯", 2) + para(doc, "扩繁批次(bre_propagation)通过两个外键形成完整链路:") + para(doc, "· scion_source_id —— 接穗来源(germplasm 或 tree 原株)", style="List Bullet") + para(doc, "· produced_clone_id —— 产出的克隆(4 棵 ramet)", style="List Bullet") + para(doc, "· 编号自描述:tree_no 形如“{品种码}-{序号}R{棵序}”,仅凭编号即可追溯血缘与棵序。", style="List Bullet") + para(doc, "追溯链:接穗原株 → 扩繁批次 → 产出克隆 → ramet 四棵 → 各棵观测。") + para(doc, "验证结果:扩繁 16,800 条,100% 有来源与产出,逐条可回查,无孤儿记录。") + heading(doc, "5.5 审定推广后的边界", 2) + para(doc, "品种审定(released)后,观测止于推广环节;生产果园的商业化栽植不再进入观测表。当前无销售/推广数据的落点(见缺口 C)。") + + # 6 性能 + heading(doc, "6 性能保障", 1) + para(doc, "记录量 318 万,表体约 990 MB,11 个索引。实测:") + table(doc, + ["场景", "实测耗时", "结论"], + [["前端按单株分页取观测", "≈ 0.39 ms", "数据库层面不是瓶颈"], + ["全性状统计取数", "≈ 802 ms", "索引足以支撑"]], + col_widths=[6.5, 3.5, 6.5]) + para(doc, "统计引擎的规模瓶颈已随稀疏重构消除:索引化设计矩阵(无稠密 Z)、稀疏 A⁻¹(无稠密 A)、共轭梯度求解,全程 O(n) 内存。个体数 ≤ N_EXACT=6000 用精确迹二分求 λ*(逐列 CG 精确 tr(A⁻¹C²²));> N_EXACT 用随机子样本(≤3000 观测 + 祖先闭包)REMl 估计方差组分(抽样 SE h² ≈ ±0.04,h²/σ² 报告值视为近似区间),再全群 MME 于子样本 λ 下解 EBV。36,720 棵树单性状 ABLUP 在 14 GiB 本机完整跑通不 OOM。") + + # 7 性状字典现状与断裂 + heading(doc, "7 性状字典现状与断裂", 1) + heading(doc, "7.1 现状", 2) + para(doc, "库里现存 20 条统计性状(分类为产量/品质/果形/营养/抗性/物候),由 ensure_reference 幂等回填;原始 43 条官方描述字典(bre_trait_seed.sql,六大类:基本鉴定/果实外观/果皮/果实大小/果肉/果核)在 wipe 重建时被清空,未回填。") + heading(doc, "7.2 断裂点", 2) + table(doc, + ["断裂点", "现象", "根因"], + [["单株评价表单", "6 个 tab 全“暂无性状”,无法录入", "表单按官方字典分类,库里只有 20 条 stats 分类"], + ["DUS", "DUS 引用 trait_code(growth_vigor 等)找不到性状", "描述字典未回填"], + ["抗病/需冷量种子", "抗病性 / 生态适应性种子性状缺失", "wipe 未回填 bre_disease_chilling_traits"], + ["通用观测下拉", "显示英文 label", "getTraitOptions 返回 label=trait_code"]], + col_widths=[3.4, 6.5, 6.5]) + para(doc, "三套来源(stats seed / 官方描述 dict / DUS-抗病 seed)互相覆盖、未统一归口,是断裂的总根因。") + + # 8 已知缺口 + heading(doc, "8 已知缺口清单", 1) + table(doc, + ["缺口", "描述", "现状", "建议方向"], + [["A 童期一年多次测定", "观测表只有 evaluate_year,一年至多一条;童期生长量真实一年多次(如春秋两次)", "无 evaluate_date", "补日期字段或次数维度"], + ["B 扩繁接穗来源手填 ID", "表单 scion_source_id 为裸数字输入,靠人工记忆原株 ID", "易错、无校验", "改为下拉,限定 germplasm/tree 来源"], + ["C 审定后无销售环节", "观测止于 released,推广量 / 销售无落点", "数据边界", "增加推广/销售记录(或明确不录)"], + ["D 通用观测下拉英文", "getTraitOptions 返回英文 trait_code 作 label", "体验差", "返回中文名"]], + col_widths=[3.2, 6.3, 3.0, 4.0]) + + # 9 待定决策 + heading(doc, "9 待定决策", 1) + table(doc, + ["决策", "选项", "说明"], + [["决策 1:标签页结构", "官方 6 类 tab / 动态 category 分组", "决定表单录入布局,需与 43 条字典匹配"], + ["决策 2:描述性状模拟值", "回填 43 条种子(含模拟值)/ 仅建字典不造值", "决定表单是否立即可用"], + ["决策 3:童期观测归属", "B 童期观测归入单株评价流程 / 独立童期页面", "决定 B 类 115 万条在哪里维护"], + ["决策 4:移动端采集入口", "本次接入 / 后续迭代", "倾向移动端采集,按评价保存"]], + col_widths=[4.2, 6.3, 6.0]) + + # 附录 + heading(doc, "附录 数据口径(2026-08-06 实测)", 1) + table(doc, + ["指标", "数值"], + [["观测总行数", "3,583,584(A 731,904 + B 1,148,160 + C 1,703,520)"], + ["A 随评构成", "45,744 个评价头 × 16 条/头 = 731,904,0 孤儿"], + ["C 克隆苗构成", "每克隆 4 棵 ramet,17 性状 × 2 年 → 136 条/克隆"], + ["扩繁批次", "16,800 条,100% 有来源与产出,无孤儿"], + ["表体规模", "约 1,129 MB / 9 索引"], + ["前端分页查询", "≈ 0.39 ms"], + ["全性状统计取数", "≈ 802 ms"], + ["模拟规模", "F1 200 组合 → 约 120k 棵树"], + ["ABLUP 全量内存", "稠密 Z(84000×86792) ≈ 54.3 GiB,苗高 ≈ 110 GiB → OOM"]], + col_widths=[6.5, 10.0]) + + doc.save(OUT) + print("OK ->", OUT) + print("size =", os.path.getsize(OUT), "bytes") + + +if __name__ == "__main__": + main() diff --git a/backend/scripts/simulate_breeding_data.py b/backend/scripts/simulate_breeding_data.py new file mode 100644 index 0000000..8147536 --- /dev/null +++ b/backend/scripts/simulate_breeding_data.py @@ -0,0 +1,1544 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +"""桃育种压测级模拟数据生成器(2026-08-06,自洽式)。 + +目标:在 dev 库(Postgres)上生成 20 年时间线的真实感育种数据—— + - ~5 万棵单株(F1 杂交圃 / F2 自交分离 / BC1 回交) + - 百万级性状观测(童期 juvenile + 成株 evaluation,含数值/文本两型) + - 完整选择晋级链(初选 sp → 复选 ap → 品系 line → 区试 regional → 审定 released), + 晋级走真实服务层 SelectionResultService(创建 clone/germplasm/pedigree、回写树阶段与状态) + - 单株评价、农事操作(采收同步产量)、定植批次、区域试验(trial/study/entry) + +设计要点: + - 参考数据(性状/人员/基地/地块/育种目标/规则/砧木/基础种质池)若缺失则补齐,不清除既有参考数据。 + - 业务数据(bre_* 表)按 FK 序全量清空重建(--no-wipe 跳过)。 + - 数值性状用简化无穷小模型:BV(父母均值+孟德尔抽样)+ 年度效应 + 误差,按性状 h² 截断到有效区间; + 选择依据 = 各性状 BV 标准化加权指数(升序性状权重取负),保证入选单株表型确实更优(选择信号真实)。 + - 过程内直调服务层(AuthSchema 超管身份,Permission 不过滤,免登录免验证码), + 批量晋级复用 create_batch(单事务逐棵 create,失败不拖垮整批)。 + - 淘汰记录(is_selected=eliminated)与淘汰树状态均在各轮正向晋级**全部完成之后**批量回写, + 避免正向晋级校验(_assert_stage_forward)先看到 eliminated 结论而误判“结论回退”。 + +运行: + ENVIRONMENT=dev python backend/scripts/simulate_breeding_data.py --dry-run # 小规模验证 + ENVIRONMENT=dev python backend/scripts/simulate_breeding_data.py --verify # 全量 + 自检 +""" +from __future__ import annotations + +import argparse +import asyncio +import json +import random +import sys +from datetime import UTC, date, datetime +from pathlib import Path + +import numpy as np + +BACKEND_DIR = Path(__file__).resolve().parents[1] +sys.path.insert(0, str(BACKEND_DIR)) + +# Windows 控制台 GBK 下中文 print 会乱码,统一以 UTF-8 输出 +try: + sys.stdout.reconfigure(encoding="utf-8", line_buffering=True) + sys.stderr.reconfigure(encoding="utf-8", line_buffering=True) +except (AttributeError, ValueError): + pass + +from sqlalchemy import func, insert, select, text # noqa: E402 +from sqlalchemy.dialects.postgresql import insert as pg_insert # noqa: E402 +from sqlalchemy.ext.asyncio import async_sessionmaker, create_async_engine # noqa: E402 + +from app.config.setting import settings # noqa: E402 +from app.core.base_schema import AuthSchema, CoreUserSchema # noqa: E402 +from app.core.bre_audit_ctx import bre_audit_suppress # noqa: E402 +from app.utils.common_util import uuid4_str # noqa: E402 + +from app.api.v1.module_bre.cross_combination.model import CrossCombinationModel # noqa: E402 +from app.api.v1.module_bre.tree.model import TreeModel # noqa: E402 +from app.api.v1.module_bre.planting.model import PlantingModel # noqa: E402 +from app.api.v1.module_bre.tree_evaluation.model import TreeEvaluationModel # noqa: E402 +from app.api.v1.module_bre.trait_observation.model import TraitObservationModel # noqa: E402 +from app.api.v1.module_bre.selection_result.model import SelectionResultModel # noqa: E402 +from app.api.v1.module_bre.field_operation.model import FieldOperationModel # noqa: E402 +from app.api.v1.module_bre.trial.model import TrialModel # noqa: E402 +from app.api.v1.module_bre.trial_study.model import TrialStudyModel # noqa: E402 +from app.api.v1.module_bre.trial_study_entry.model import TrialStudyEntryModel # noqa: E402 +from app.api.v1.module_bre.trait.model import TraitModel # noqa: E402 +from app.api.v1.module_bre.germplasm.model import BreedingGermplasmModel # noqa: E402 +from app.api.v1.module_bre.target.model import TargetModel # noqa: E402 +from app.api.v1.module_bre.plan.model import PlanModel # noqa: E402 +from app.api.v1.module_bre.personnel.model import PersonnelModel # noqa: E402 +from app.api.v1.module_bre.site.model import BreedingSiteModel, BreedingPlotModel # noqa: E402 +from app.api.v1.module_bre.rootstock.model import RootstockModel # noqa: E402 +from app.api.v1.module_bre.selection_rule.model import SelectionRuleModel # noqa: E402 +from app.api.v1.module_bre.selection_result.service import SelectionResultService # noqa: E402 +from app.api.v1.module_bre.selection_result.schema import SelectionResultBatchCreateSchema # noqa: E402 +from app.api.v1.module_bre.pollination.model import PollinationModel # noqa: E402 +from app.api.v1.module_bre.pollen.model import PollenModel # noqa: E402 +from app.api.v1.module_bre.seed_lot.model import SeedLotModel # noqa: E402 +from app.api.v1.module_bre.seed_treatment.model import SeedTreatmentModel # noqa: E402 +from app.api.v1.module_bre.seedling.model import SeedlingModel # noqa: E402 +from app.api.v1.module_bre.propagation.model import PropagationModel # noqa: E402 +from app.api.v1.module_bre.environment_condition.model import EnvironmentConditionModel # noqa: E402 +from app.api.v1.module_bre.observation.model import ObservationModel # noqa: E402 +from app.api.v1.module_bre.treatment.model import TreatmentModel # noqa: E402 +from app.api.v1.module_bre.clone.model import CloneModel # noqa: E402 +# UserModel 必须显式导入:多个 module_bre 模型的 relationship 以字符串引用 UserModel, +# 未注册进 mapper registry 会在首次配置映射时抛 InvalidRequestError +from app.api.v1.module_system.user.model import UserModel # noqa: E402,F401 + +# --------------------------------------------------------------------------- +# 常量与参数 +# --------------------------------------------------------------------------- +WINDOW_END = 2026 # 数据“今天”年份(2026-08-06) +F1_YEARS = range(2007, 2024) # F1 杂交组合年份 2007..2023 +F1_TREES = 60 # 10 人团队真实规模:F1 = 17年×12×60 = 12,240 +F2_COMBOS = 120 +F2_TREES = 30 # F2 = 120×30 = 3,600 +BC_COMBOS = 40 +BC_TREES = 20 # BC = 40×20 = 800(全圃活体合计 ≈ 3.4 万,含 ramet) +K_RAMETS = 4 # 每入选克隆的无性系苗(ramet)数:保证 ABLUP G1 每 clone ≥2 样本 +RAMET_EVAL_YEARS = 2 # 无性系苗成株/童期观测年数(自初选年起的连续年份) + +# 晋级轮次:名称 / 选育结论 / 入选比例 / 起始阶段 / 晋级阶段 / 距定植年偏移 +ROUNDS = [ + ("sp", "selected", 0.35, "seedling", "sp", 3), + ("ap", "primary", 0.50, "sp", "ap", 5), + ("line", "key", 0.40, "ap", "line", 6), + ("regional", "preserved", 0.50, "line", "regional_trial", 7), + ("released", "preserved", 0.30, "regional_trial", "released", 8), +] + +# 数值性状:code, name, unit, category, vmin, vmax, stage, direction, into_ebv, mu, sd, h2, weight +NUMERIC_TRAITS = [ + ("single_tree_yield", "单株产量", "kg", "产量", 0, 60, "evaluation", "desc", "1", 18.0, 8.0, 0.35, 0.25), + ("fruit_number", "结果个数", "个", "产量", 0, 600, "evaluation", "desc", "1", 220.0, 80.0, 0.30, 0.03), + ("avg_fruit_weight", "平均单果重", "g", "果形", 50, 400, "evaluation", "desc", "1", 210.0, 50.0, 0.40, 0.20), + ("marketable_rate", "好果率", "%", "品质", 0, 100, "evaluation", "desc", "1", 82.0, 12.0, 0.30, 0.15), + ("ssc", "可溶性固形物", "%", "品质", 5, 22, "evaluation", "desc", "1", 12.5, 2.2, 0.55, 0.20), + ("firmness_value", "果肉硬度", "kg/cm²", "品质", 0, 20, "evaluation", "desc", "1", 8.5, 2.5, 0.40, 0.05), + ("acidity", "可滴定酸", "%", "品质", 0, 3.5, "evaluation", "asc", "1", 0.7, 0.25, 0.40, -0.05), + ("red_blush", "着色度", "%", "果形", 0, 100, "evaluation", "desc", "1", 55.0, 22.0, 0.50, 0.05), + ("chilling_requirement", "需冷量", "h", "物候", 0, 1200, "evaluation", "asc", "1", 750.0, 180.0, 0.60, -0.05), + ("stone_weight", "核重", "g", "果形", 0, 30, "evaluation", "asc", "1", 8.0, 3.0, 0.30, -0.02), + ("fruit_length", "果实纵径", "mm", "果形", 40, 110, "evaluation", "desc", "1", 68.0, 10.0, 0.35, 0.02), + ("fruit_diameter", "果实横径", "mm", "果形", 40, 120, "evaluation", "desc", "1", 72.0, 12.0, 0.35, 0.02), + ("uniformity", "果实整齐度", "%", "品质", 0, 100, "evaluation", "desc", "1", 78.0, 10.0, 0.25, 0.02), + ("seedling_height", "苗高", "cm", "营养", 0, 400, "juvenile", "desc", "0", 180.0, 50.0, 0.30, 0.0), + ("trunk_diameter", "地径", "cm", "营养", 0, 20, "juvenile", "desc", "0", 6.0, 1.8, 0.30, 0.0), + ("shoot_growth", "新梢生长量", "cm", "营养", 0, 250, "juvenile", "desc", "0", 90.0, 35.0, 0.25, 0.0), + ("leaf_disease_severity", "叶部病害程度", "%", "抗性", 0, 100, "juvenile", "asc", "0", 15.0, 12.0, 0.20, 0.0), +] + +# 文本性状:code, name, category, stage, options, probs +TEXT_TRAITS = [ + ("disease_resistance", "抗病性", "抗性", "evaluation", ["高抗", "抗", "中抗", "感"], [0.12, 0.35, 0.35, 0.18]), + ("flesh_color", "果肉颜色", "品质", "evaluation", ["白色", "黄白", "黄色", "红色"], [0.35, 0.25, 0.30, 0.10]), + ("fruit_shape", "果实形状", "果形", "evaluation", ["圆形", "扁圆", "蟠桃形", "椭圆形"], [0.40, 0.20, 0.15, 0.25]), +] + +# 描述性状(单株鉴定表单 6 标签页官方字典,权威来源:育种管理系统操作手册)。 +# 结构:(code, name, category, data_type, options_or_None, mu_or_None, sd_or_None, vmin_or_None, vmax_or_None, reuse_stats_code_or_None) +# data_type: categorical(分级,选项入 scale_json) / date / numeric +# reuse_stats_code: 值同源于某 stats 观测(如纵径横径复用 stats fruit_length/fruit_diameter),不独立造值 +# 归并策略:avg_fruit_weight/ssc/flesh_color/fruit_shape 与 stats 撞 code → 描述侧定义覆盖升级(category/data_type/scale_json); +# 纵径/横径(longitudinal_dia/transverse_dia)与 stats fruit_length/fruit_diameter 语义重复 → 造值时取同源值。 +DESCRIPTIVE_TRAITS = [ + # ① 基本鉴定 + ("growth_vigor", "生长势", "基本鉴定", "categorical", ["强", "中", "弱"], None, None, None, None, None), + ("growth_type", "生长型", "基本鉴定", "categorical", ["普通", "半矮", "矮"], None, None, None, None, None), + ("first_flower_date", "始花日期", "基本鉴定", "date", None, None, None, None, None, None), + ("flower_type", "花型", "基本鉴定", "categorical", ["蔷薇型", "铃型"], None, None, None, None, None), + ("pollen", "花粉", "基本鉴定", "categorical", ["无", "少", "中", "多"], None, None, None, None, None), + ("load_amount", "负载量", "基本鉴定", "categorical", ["高", "中", "低"], None, None, None, None, None), + ("maturity_uniformity", "成熟一致性", "基本鉴定", "categorical", ["一致", "较一致", "不一致"], None, None, None, None, None), + # ② 果实外观 + ("fruit_shape", "果形", "果实外观", "categorical", ["扁平", "扁圆", "圆", "椭圆", "卵圆"], None, None, None, None, None), + ("apex_shape", "果顶", "果实外观", "categorical", ["突出", "稍突出", "圆平", "稍凹陷", "凹陷"], None, None, None, None, None), + ("stem_cavity_depth", "梗洼深度", "果实外观", "categorical", ["浅", "中", "深"], None, None, None, None, None), + ("stem_cavity_width", "梗洼宽度", "果实外观", "categorical", ["窄", "中", "宽"], None, None, None, None, None), + ("fruit_base", "果基", "果实外观", "categorical", ["正", "稍偏", "偏"], None, None, None, None, None), + ("suture_line", "缝合线", "果实外观", "categorical", ["浅", "中", "深"], None, None, None, None, None), + ("symmetry", "对称性", "果实外观", "categorical", ["对称", "较对称", "不对称"], None, None, None, None, None), + # ③ 果皮 + ("pubescence", "茸毛", "果皮", "categorical", ["无", "少", "中", "多"], None, None, None, None, None), + ("skin_base_color", "底色", "果皮", "categorical", ["淡绿", "绿白", "白", "浅黄", "深黄"], None, None, None, None, None), + ("coloration_area", "着色面积", "果皮", "categorical", ["无", "少", "中", "多", "全"], None, None, None, None, None), + ("coloration_depth", "着色深度", "果皮", "categorical", ["浅", "中", "深"], None, None, None, None, None), + ("coloration_brightness", "着色亮度", "果皮", "categorical", ["暗", "中", "亮"], None, None, None, None, None), + ("coloration_pattern", "着色形态", "果皮", "categorical", ["晕", "条", "斑"], None, None, None, None, None), + ("skin_spot", "果面斑点", "果皮", "categorical", ["少", "中", "多"], None, None, None, None, None), + ("skin_peelability", "果皮剥离度", "果皮", "categorical", ["不能", "难", "易"], None, None, None, None, None), + # ④ 果实大小(数值 core) + ("max_fruit_weight", "最大果重", "果实大小", "numeric", None, 300.0, 70.0, 100, 600, None), + ("avg_fruit_weight", "平均果重", "果实大小", "numeric", None, 210.0, 50.0, 50, 400, None), + ("longitudinal_dia", "纵径", "果实大小", "numeric", None, None, None, None, None, "fruit_length"), + ("transverse_dia", "横径", "果实大小", "numeric", None, None, None, None, None, "fruit_diameter"), + ("lateral_dia", "侧径", "果实大小", "numeric", None, 70.0, 11.0, 40, 115, None), + # ⑤ 果肉 + ("flesh_thickness", "果肉厚度", "果肉", "numeric", None, 18.0, 3.0, 5, 35, None), + ("flesh_color", "肉色", "果肉", "categorical", ["淡绿", "白", "黄白", "黄", "橙黄", "红", "紫红"], None, None, None, None, None), + ("flesh_firmness", "硬度", "果肉", "categorical", ["很软", "软", "中", "硬", "很硬"], None, None, None, None, None), + ("subepidermal_anthocyanin", "皮下花色苷", "果肉", "categorical", ["无", "少", "中", "多"], None, None, None, None, None), + ("flesh_anthocyanin", "果肉花色苷", "果肉", "categorical", ["无", "少", "中", "多"], None, None, None, None, None), + ("stone_cavity_anthocyanin", "核窝花色苷", "果肉", "categorical", ["无", "少", "中", "多"], None, None, None, None, None), + ("ssc", "SSC", "果肉", "numeric", None, 12.5, 2.2, 5, 22, None), + ("fiber", "纤维", "果肉", "categorical", ["少", "中", "多"], None, None, None, None, None), + ("juice", "汁液", "果肉", "categorical", ["少", "中", "多"], None, None, None, None, None), + ("sweet_acidity", "甜酸度", "果肉", "categorical", ["淡甜", "甜", "浓甜", "酸甜", "甜酸适中", "甜酸", "酸"], None, None, None, None, None), + ("aroma", "香气", "果肉", "categorical", ["淡", "中", "多"], None, None, None, None, None), + ("quality_overall", "品质综合", "果肉", "categorical", ["上", "中", "下"], None, None, None, None, None), + # ⑥ 果核 + ("stone_adherence", "粘离性", "果核", "categorical", ["粘核", "半离", "离核"], None, None, None, None, None), + ("stone_size", "核大小", "果核", "categorical", ["小", "中", "大"], None, None, None, None, None), + ("stone_shape", "核形", "果核", "categorical", ["扁平", "近圆", "椭圆", "倒卵圆", "卵圆"], None, None, None, None, None), + ("stone_cracking", "裂核", "果核", "categorical", ["无", "少", "中", "多"], None, None, None, None, None), +] +# 描述性状中与 stats 撞 code、由描述侧覆盖升级现有行的 code 集 +DESCRIPTIVE_OVERRIDE_CODES = {"avg_fruit_weight", "ssc", "flesh_color", "fruit_shape"} +# 描述性状中值同源复用 stats 观测的映射 {描述code: stats code} +DESCRIPTIVE_REUSE_STATS = {t[0]: t[9] for t in DESCRIPTIVE_TRAITS if t[9]} +# 描述性状 numeric core 的单位与先验遗传力(抄 bre_trait_seed.sql 官方值) +DESCRIPTIVE_UNITS = {"max_fruit_weight": "g", "avg_fruit_weight": "g", "longitudinal_dia": "mm", + "transverse_dia": "mm", "lateral_dia": "mm", "flesh_thickness": "mm", "ssc": "%"} +DESCRIPTIVE_H2 = {"max_fruit_weight": 0.5, "avg_fruit_weight": 0.5, "longitudinal_dia": 0.4, + "transverse_dia": 0.4, "lateral_dia": 0.4, "flesh_thickness": 0.3, "ssc": 0.5} +# 描述 numeric core 的遗传代理(不在 bv 矩阵列中):{描述code: (stats code, 系数)}, +# 借同源 stats 性状的 BV 保持遗传相关(最大果重↔平均果重、侧径/肉厚↔横径) +DESCRIPTIVE_BV_PROXY = { + "max_fruit_weight": ("avg_fruit_weight", 1.3), + "lateral_dia": ("fruit_diameter", 0.9), + "flesh_thickness": ("fruit_diameter", 0.3), +} + +# 基础种质池(F1 亲本):name, variety_type, maturity, firmness, s_alleles, avg_fruit_weight, ssc, chilling +FOUNDERS = [ + ("春美", "普通桃", "早", "硬溶质", "S1/S2", 150, 11.5, 650), + ("春雪", "普通桃", "早", "硬溶质", "S1/S3", 180, 12.0, 700), + ("白凤", "普通桃", "中", "软溶质", "S2/S3", 220, 13.0, 780), + ("湖景蜜露", "普通桃", "中", "软溶质", "S3/S4", 240, 12.5, 800), + ("京玉", "普通桃", "中", "硬溶质", "S1/S4", 260, 12.0, 750), + ("大久保", "普通桃", "中", "软溶质", "S2/S4", 280, 11.0, 820), + ("仓方早生", "普通桃", "早", "硬溶质", "S1/S9", 190, 10.5, 700), + ("冈山白", "普通桃", "中", "软溶质", "S2/S9", 300, 12.0, 850), + ("早香玉", "普通桃", "早", "硬溶质", "S3/S9", 160, 11.0, 680), + ("雨花露", "普通桃", "早", "软溶质", "S1/S7", 170, 10.0, 690), + ("霞晖1号", "普通桃", "早", "硬溶质", "S2/S7", 175, 11.5, 670), + ("砂子早生", "普通桃", "早", "硬溶质", "S3/S7", 165, 10.8, 660), + ("中油桃4号", "油桃", "中", "硬溶质", "S1/S3", 210, 12.8, 760), + ("中油桃8号", "油桃", "中", "硬溶质", "S2/S4", 235, 13.5, 780), + ("瑞光18", "油桃", "中", "硬溶质", "S1/S4", 220, 12.2, 770), + ("瑞光28", "油桃", "晚", "硬溶质", "S3/S9", 245, 13.0, 820), + ("瑞蟠4号", "蟠桃", "中", "软溶质", "S2/S3", 210, 12.5, 800), + ("瑞蟠21", "蟠桃", "中", "硬溶质", "S1/S2", 230, 12.8, 790), + ("晚白桃", "普通桃", "晚", "硬溶质", "S7/S9", 270, 12.5, 900), + ("京红", "普通桃", "中", "软溶质", "S1/S9", 250, 11.8, 790), + ("燕红", "普通桃", "中", "硬溶质", "S2/S9", 265, 12.2, 810), + ("新川中岛", "普通桃", "中", "硬溶质", "S3/S4", 280, 13.2, 830), + ("川中岛白桃", "普通桃", "晚", "硬溶质", "S1/S7", 290, 13.5, 860), + ("锦绣黄桃", "黄肉", "中", "硬溶质", "S2/S7", 260, 14.0, 800), + ("锦园", "黄肉", "中", "硬溶质", "S1/S4", 275, 13.8, 820), + ("金童5号", "黄肉", "中", "硬溶质", "S3/S7", 240, 13.0, 790), + ("金童8号", "黄肉", "中", "硬溶质", "S2/S9", 250, 13.2, 800), + ("早黄金", "黄肉", "早", "硬溶质", "S1/S3", 200, 12.8, 720), + ("中蟠桃1号", "蟠桃", "中", "硬溶质", "S2/S4", 225, 12.0, 790), + ("中蟠桃11", "蟠桃", "晚", "硬溶质", "S1/S9", 235, 12.6, 830), + ("中油蟠1号", "油蟠桃", "中", "硬溶质", "S3/S4", 230, 12.8, 810), + ("中华寿桃", "普通桃", "晚", "硬溶质", "S2/S9", 320, 13.0, 950), + ("大五星", "普通桃", "中", "硬溶质", "S1/S4", 300, 12.2, 840), + ("布目早生", "普通桃", "早", "软溶质", "S3/S9", 195, 10.8, 710), + ("冬桃", "普通桃", "晚", "硬溶质", "S7/S9", 310, 13.5, 1100), + ("血桃", "普通桃", "中", "硬溶质", "S1/S7", 200, 11.0, 780), + ("脆桃", "普通桃", "早", "硬溶质", "S2/S3", 175, 10.5, 690), + ("黄肉晚熟", "黄肉", "晚", "硬溶质", "S1/S2", 280, 14.2, 950), + ("白桃(地方)", "普通桃", "中", "软溶质", "S2/S7", 260, 11.5, 850), + ("油桃晚熟", "油桃", "晚", "硬溶质", "S3/S7", 255, 13.8, 880), +] + +PERSONNEL = [ + ("王育桃", "男", "育种家", "主持杂交组合与选择晋级"), + ("李农艺", "女", "农艺师", "负责定植与田间评价"), + ("张田间", "男", "技术员", "负责农事操作与采收记录"), +] +SITES = [("主育种基地", "河南郑州", "华北", 35.6, 113.6, 150, "核心育种圃与品比圃"), + ("试验示范站", "山东泰安", "华北", 36.2, 117.1, 250, "区试与示范栽培点")] +PLOTS = [(0, "P01", "南北行", "选种圃"), + (0, "P02", "南北行", "选种圃"), + (0, "P03", "南北行", "复选圃"), + (0, "P04", "东西行", "品比圃"), + (0, "P05", "东西行", "育种圃"), + (0, "P06", "南北行", "杂种圃"), + (1, "T01", "南北行", "区试圃"), + (1, "T02", "东西行", "区试圃"), + (1, "T03", "南北行", "区试圃")] +PLANS = [("桃-2020-1", "桃树品种选育计划(2020-2030)", "早熟优质大果硬肉耐贮综合育种")] +TARGETS = [("早熟优质", "早熟", "早熟、优质、丰产,树势开张"), + ("大果硬肉", "硬肉", "大果型、硬溶质、耐贮运"), + ("抗病耐贮", "抗病", "抗细菌性穿孔病、耐贮运"), + ("低需冷量", "冷量", "低需冷量、适南岭以南促早栽培")] +RULES = [("初选", "sp", "综合评价得分与指数筛选"), + ("复选", "ap", "连续两年结果表现复核"), + ("品系比较", "line", "品系比较试验表现"), + ("区域试验", "regional_trial", "区域试验产量与品质表现")] +ROOTSTOCKS = [("毛桃", "乔化", "强"), ("山桃", "乔化", "中"), ("GF677", "半矮化", "强")] + + +def _now() -> datetime: + return datetime.now(UTC) + + +def _base() -> dict: + return {"uuid": uuid4_str(), "created_time": _now(), "updated_time": _now(), "is_deleted": False} + + +# --------------------------------------------------------------------------- +# 工具 +# --------------------------------------------------------------------------- +def chunked(seq: list, size: int): + for i in range(0, len(seq), size): + yield seq[i:i + size] + + +def clip_round(v, vmin: float, vmax: float): + return float(np.clip(np.round(v, 1), vmin, vmax)) + + +async def insert_chunked(session, model, rows: list[dict], size: int = 20000) -> int: + total = 0 + for chunk in chunked(rows, size): + await session.execute(insert(model), chunk) + total += len(chunk) + return total + + +# --------------------------------------------------------------------------- +# 参考数据补齐 +# --------------------------------------------------------------------------- +# 官方字典种子(产量/品质/物候/生长量/描述 + 抗病/需冷量 + DUS TG/53 描述符)以 +# backend/sql/*.sql 为单一事实源。模拟脚本只为其自有性状造观测,不建这些字典行; +# wipe 重建若不同步回填,回归断言(e2e_busflow 19 码 / e2e_refine DUS)必缺。 +# 幂等:各文件均为单个 INSERT ... ON CONFLICT(DO UPDATE / DO NOTHING)。 +_REFERENCE_SQL_SEEDS = [ + "bre_yield_phenology_traits.sql", + "bre_disease_chilling_traits.sql", + "bre_dus_seed.sql", +] + + +async def _apply_reference_sql_seeds(session) -> None: + sql_dir = Path(__file__).resolve().parent.parent / "sql" + for fname in _REFERENCE_SQL_SEEDS: + content = (sql_dir / fname).read_text(encoding="utf-8") + for stmt in content.split(";"): + body = "\n".join(ln for ln in stmt.splitlines() + if not ln.strip().startswith("--")).strip() + if body: + await session.execute(text(body)) + await session.commit() + + +async def ensure_reference(session) -> dict: + """补齐参考数据(缺失才插入),返回供后续使用的 id 映射。""" + ref = {} + + # 性状(trait_code 幂等) + existing_codes = {r[0] for r in await session.execute(select(TraitModel.trait_code))} + trait_upserts = [] + for code, name, unit, cat, vmin, vmax, stage, direction, into_ebv, *_ in NUMERIC_TRAITS: + if code not in existing_codes: + trait_upserts.append({"trait_code": code, "trait_name": name, "category": cat, "data_type": "numeric", + "unit": unit, "is_core": "1", "valid_min": vmin, "valid_max": vmax, + "stage": stage, "direction": direction, "into_ebv": into_ebv, + "default_h2": next(t[11] for t in NUMERIC_TRAITS if t[0] == code)}) + for code, name, cat, stage, options, probs in TEXT_TRAITS: + if code not in existing_codes: + trait_upserts.append({"trait_code": code, "trait_name": name, "category": cat, "data_type": "text", + "unit": None, "is_core": "1", "valid_min": None, "valid_max": None, + "stage": stage, "direction": "desc", "into_ebv": "0", + "scale_json": json.dumps(list(zip(options, probs)), ensure_ascii=False), + "default_h2": 0.2}) + if trait_upserts: + await session.execute(insert(TraitModel), [{**_base(), **r} for r in trait_upserts]) + await session.commit() + # 描述性状(43 条官方字典):撞 code 的覆盖升级(category/data_type/scale_json 对齐描述侧), + # 其余幂等插入。desc 覆盖不建重复行(avg_fruit_weight/ssc/flesh_color/fruit_shape 与 stats 同归口)。 + desc_upserts = [] + for code, name, cat, dtype, options, mu, sd, vmin, vmax, _reuse in DESCRIPTIVE_TRAITS: + if dtype == "numeric": + desc_upserts.append({"trait_code": code, "trait_name": name, "category": cat, "data_type": "numeric", + "unit": DESCRIPTIVE_UNITS.get(code), "is_core": "1", + "valid_min": vmin, "valid_max": vmax, "stage": "evaluation", + "direction": "desc", "into_ebv": "1", + "scale_json": None, "default_h2": DESCRIPTIVE_H2.get(code)}) + elif dtype == "categorical": + desc_upserts.append({"trait_code": code, "trait_name": name, "category": cat, "data_type": "categorical", + "unit": None, "is_core": "0", "valid_min": None, "valid_max": None, + "stage": "evaluation", "direction": "desc", "into_ebv": "0", + "scale_json": json.dumps(options, ensure_ascii=False), "default_h2": None}) + else: # date + desc_upserts.append({"trait_code": code, "trait_name": name, "category": cat, "data_type": "date", + "unit": None, "is_core": "0", "valid_min": None, "valid_max": None, + "stage": "evaluation", "direction": "desc", "into_ebv": "0", + "scale_json": None, "default_h2": None}) + for r in desc_upserts: + await session.execute( + pg_insert(TraitModel) + .values({**_base(), **r}) + .on_conflict_do_update( + index_elements=["trait_code"], + # is_core 不入 set_:撞 code 的 stats 行(flesh_color/fruit_shape 等)保留原 is_core, + # 仅对齐 category/data_type/scale_json 等描述侧定义 + set_={"trait_name": r["trait_name"], "category": r["category"], "data_type": r["data_type"], + "unit": r["unit"], "valid_min": r["valid_min"], + "valid_max": r["valid_max"], "scale_json": r["scale_json"], + "direction": r["direction"], "into_ebv": r["into_ebv"], "default_h2": r["default_h2"], + "updated_time": func.now()})) + await session.commit() + await _apply_reference_sql_seeds(session) + rows = await session.execute(select(TraitModel.id, TraitModel.trait_code)) + ref["trait_id_by_code"] = {code: tid for tid, code in rows.all()} + + # 人员 + if (await session.scalar(select(func.count(PersonnelModel.id))) or 0) == 0: + await session.execute(insert(PersonnelModel), [ + {**_base(), "name": n, "gender": g, "role": r, "remark": rm} + for n, g, r, rm in PERSONNEL + ]) + await session.commit() + ref["personnel_ids"] = [r[0] for r in await session.execute( + select(PersonnelModel.id).where(PersonnelModel.is_deleted == False))] # noqa: E712 + + # 基地/地块 + if (await session.scalar(select(func.count(BreedingSiteModel.id))) or 0) == 0: + await session.execute(insert(BreedingSiteModel), [ + {**_base(), "site_name": n, "address": a, "eco_type": e, "latitude": la, "longitude": lo, + "elevation": el, "remark": rm} + for n, a, e, la, lo, el, rm in SITES + ]) + await session.commit() + site_ids = [r[0] for r in await session.execute(select(BreedingSiteModel.id).order_by(BreedingSiteModel.id))] + if (await session.scalar(select(func.count(BreedingPlotModel.id))) or 0) == 0: + await session.execute(insert(BreedingPlotModel), [ + {**_base(), "site_id": site_ids[si], "plot_code": code, "row_orientation": ro, "grid_note": note, + "row_count": 40, "col_count": 30, "area": 5.0} + for si, code, ro, note in PLOTS + ]) + await session.commit() + plot_rows = await session.execute( + select(BreedingPlotModel.id, BreedingPlotModel.site_id).order_by(BreedingPlotModel.id) + .where(BreedingPlotModel.is_deleted == False)) # noqa: E712 + ref["plots"] = [{"id": pid, "site_id": sid} for pid, sid in plot_rows.all()] + + # 计划/目标/规则 + if (await session.scalar(select(func.count(PlanModel.id))) or 0) == 0: + await session.execute(insert(PlanModel), [ + {**_base(), "plan_code": c, "plan_name": n, "objective": o} + for c, n, o in PLANS + ]) + await session.commit() + plan_ids = [r[0] for r in await session.execute( + select(PlanModel.id).where(PlanModel.is_deleted == False))] # noqa: E712 + ref["plan_ids"] = plan_ids + if (await session.scalar(select(func.count(TargetModel.id))) or 0) == 0: + await session.execute(insert(TargetModel), [ + {**_base(), "target_name": n, "series": s, "description": d, "is_preset": "1", + "plan_id": plan_ids[0] if plan_ids else None} + for n, s, d in TARGETS + ]) + await session.commit() + target_rows = await session.execute( + select(TargetModel.id).order_by(TargetModel.id).where(TargetModel.is_deleted == False)) # noqa: E712 + ref["target_ids"] = [r[0] for r in target_rows.all()] + if (await session.scalar(select(func.count(SelectionRuleModel.id))) or 0) == 0: + target0 = ref["target_ids"][0] if ref["target_ids"] else None + await session.execute(insert(SelectionRuleModel), [ + {**_base(), "rule_name": n, "stage": st, "target_id": target0, "conditions_json": "[]", + "logic": "and", "action": "promote", "priority": 1, "enabled": "1", "remark": rm} + for n, st, rm in RULES + ]) + await session.commit() + rule_rows = await session.execute( + select(SelectionRuleModel.id).order_by(SelectionRuleModel.id) + .where(SelectionRuleModel.is_deleted == False)) # noqa: E712 + ref["rule_ids"] = [r[0] for r in rule_rows.all()] + + # 砧木 + if (await session.scalar(select(func.count(RootstockModel.id))) or 0) == 0: + await session.execute(insert(RootstockModel), [ + {**_base(), "rootstock_name": n, "dwarf_class": d, "compatibility": c} + for n, d, c in ROOTSTOCKS + ]) + await session.commit() + root_rows = await session.execute( + select(RootstockModel.id).order_by(RootstockModel.id) + .where(RootstockModel.is_deleted == False)) # noqa: E712 + ref["rootstock_ids"] = [r[0] for r in root_rows.all()] + + # 基础种质池(按名称幂等补齐) + germ_rows = await session.execute(select(BreedingGermplasmModel.id, BreedingGermplasmModel.cultivar_name)) + germ_by_name = {name: gid for gid, name in germ_rows.all()} + founder_missing = [FOUNDERS[i] for i, (name, *_rest) in enumerate(FOUNDERS) if name not in germ_by_name] + if founder_missing: + await session.execute(insert(BreedingGermplasmModel), [ + {**_base(), "cultivar_name": n, "variety_type": vt, "maturity_period": m, "firmness": f, + "s_alleles": s, "avg_fruit_weight": w, "ssc": ssc, "chilling_requirement": c, + "can_be_female": True, "can_be_male": True, "stage": "germplasm", "generation": "F0", + "storage_type": "田间", "is_rootstock": False, "biological_status": "landrace", + "breeding_program": "桃育种项目"} + for n, vt, m, f, s, w, ssc, c in founder_missing + ]) + await session.commit() + germ_rows = await session.execute(select(BreedingGermplasmModel.id, BreedingGermplasmModel.cultivar_name)) + germ_by_name = {name: gid for gid, name in germ_rows.all()} + ref["founder_ids"] = [germ_by_name[f[0]] for f in FOUNDERS] + return ref + + +# --------------------------------------------------------------------------- +# 清表 +# --------------------------------------------------------------------------- +async def wipe_business(session) -> None: + """清空全部 bre_* 业务表(含参考数据与历史测试残留),保留字典(sys_dict_*)。 + + 动态取当前库中所有 bre_ 前缀表,TRUNCATE ... RESTART IDENTITY CASCADE 一次清掉 + 整个外键依赖闭包(已确认无非 bre_ 表外键指向 bre_ 表),确保最终库中 + 只保留字典 + 本脚本重新生成的数据。 + """ + tbls = [r[0] for r in await session.execute(text( + "SELECT tablename FROM pg_tables WHERE schemaname='public' AND tablename LIKE 'bre_%'"))] + if tbls: + await session.execute(text(f"TRUNCATE TABLE {', '.join(tbls)} RESTART IDENTITY CASCADE")) + await session.commit() + + +# --------------------------------------------------------------------------- +# 遗传模拟与选择 +# --------------------------------------------------------------------------- +def build_gen_params() -> tuple[list[str], dict, np.ndarray, np.ndarray]: + """返回 (trait_codes, param_map, sigmaA, w)。""" + codes = [t[0] for t in NUMERIC_TRAITS] + params = {} + sigmaA = np.zeros(len(codes)) + w = np.zeros(len(codes)) + for i, t in enumerate(NUMERIC_TRAITS): + code, _name, _unit, _cat, vmin, vmax, _stage, direction, into_ebv, mu, sd, h2, weight = t + params[code] = {"mu": mu, "sd": sd, "h2": h2, "vmin": vmin, "vmax": vmax, + "direction": direction, "into_ebv": into_ebv, "stage": _stage, "idx": i} + sigmaA[i] = sd * np.sqrt(h2) + w[i] = weight + wsum = np.abs(w).sum() + if wsum > 0: + w = w / wsum + return codes, params, sigmaA, w + + +def idx_score(bv: np.ndarray, sd: np.ndarray, w: np.ndarray) -> np.ndarray: + """加权标准化 BV 指数(高 = 优)。""" + return (bv / sd) @ w + + +def sample_bv_parents(rng: np.random.Generator, dam_bv: np.ndarray, sire_bv: np.ndarray, sigmaA: np.ndarray) -> np.ndarray: + """后代 BV = 中亲 + 孟德尔抽样(sd = sigmaA/sqrt(2))。""" + return 0.5 * (dam_bv + sire_bv) + rng.normal(0.0, 1.0, size=sigmaA.shape) * (sigmaA * 0.7071) + + +def select_pool(pool: np.ndarray, idx: np.ndarray, frac: float, rng: np.random.Generator) -> tuple[np.ndarray, np.ndarray]: + """在候选池内按指数取前 frac(含小随机扰动保持多样性),返回 (选中, 淘汰)。""" + n = len(pool) + k = max(1, int(round(n * frac))) + if k >= n: + return pool, np.array([], dtype=int) + jitter = rng.normal(0.0, 0.05, size=n) + rank = np.argsort(-(idx[pool] + jitter)) + return pool[rank[:k]], pool[rank[k:]] + + +async def apply_elim_status(session) -> None: + """正向晋级全部完成后,将存在淘汰记录的树置为终止态(只进不退,安全)。""" + await session.execute(text( + "UPDATE bre_tree t SET status='eliminated' " + "WHERE EXISTS (SELECT 1 FROM bre_selection_result s " + "WHERE s.tree_id = t.id AND s.is_deleted = false AND s.is_selected = 'eliminated')")) + + +# --------------------------------------------------------------------------- +# 主流程 +# --------------------------------------------------------------------------- +async def main() -> None: + parser = argparse.ArgumentParser(description="桃育种压测级模拟数据生成器") + parser.add_argument("--dry-run", action="store_true", help="小规模验证(精简组合/株数)") + parser.add_argument("--seed", type=int, default=42, help="随机种子") + parser.add_argument("--no-wipe", action="store_true", help="不清空既有业务数据") + parser.add_argument("--verify", action="store_true", help="生成后运行自检") + parser.add_argument("--batch-size", type=int, default=50, help="批量晋级每批单株数") + parser.add_argument("--end-year", type=int, default=WINDOW_END, help="数据终止年份") + args = parser.parse_args() + + if settings.DATABASE_TYPE != "postgres": + raise SystemExit("本脚本仅支持 Postgres(dev 环境)。请以 ENVIRONMENT=dev 运行。") + + dry = args.dry_run + f1_years = list(F1_YEARS) if not dry else [2007, 2008, 2009, 2010] + f1_per_year = 3 if dry else 12 + f1_trees = 24 if dry else F1_TREES + f2_combos = 3 if dry else F2_COMBOS + f2_trees = 16 if dry else F2_TREES + bc_combos = 2 if dry else BC_COMBOS + bc_trees = 12 if dry else BC_TREES + end_year = args.end_year + + rng = np.random.default_rng(args.seed) + pyrand = random.Random(args.seed) + codes, params, sigmaA, w = build_gen_params() + nT = len(codes) + sd = np.array([params[c]["sd"] for c in codes]) + + engine = create_async_engine(settings.ASYNC_DB_URI, pool_size=8, max_overflow=12, pool_pre_ping=True) + Session = async_sessionmaker(engine, expire_on_commit=False) + auth = AuthSchema(user=CoreUserSchema(id=1, username="admin", name="管理员", is_superuser=True)) + t0 = datetime.now(UTC) + counts = {"trees": 0, "evals": 0, "obs": 0, "field_ops": 0, "selections": 0, + "combos": 0, "plantings": 0, "clones": 0, "germplasm": 0, "pedigree": 0, + "trials": 0, "studies": 0, "entries": 0, + "pollination": 0, "pollen": 0, "seed_lots": 0, "seed_treatments": 0, "seedlings": 0, + "propagations": 0, "env_conds": 0, "observations": 0, "treatments": 0} + print(f"[sim] ENVIRONMENT={settings.ENVIRONMENT} DB={settings.DATABASE_TYPE}@{settings.DATABASE_NAME} dry_run={dry}") + + async with Session() as session: + if not args.no_wipe: + print("[sim] 1/8 清空全部 bre_* 业务表(字典保留)...") + await wipe_business(session) + print("[sim] 2/8 补齐参考数据(性状/人员/基地/地块/目标/规则/砧木/基础种质池)...") + ref = await ensure_reference(session) + trait_id = ref["trait_id_by_code"] + personnel_ids = ref["personnel_ids"] + plots = ref["plots"] + founder_ids = ref["founder_ids"] + nf = len(founder_ids) + rule_ids = ref["rule_ids"] + rootstock_ids = ref["rootstock_ids"] + target_ids = ref["target_ids"] + plan_ids = ref["plan_ids"] + + # 基础种质池遗传值(标准正态 × 加性标准差) + founder_bv = rng.normal(0.0, 1.0, (nf, nT)) * sigmaA[None, :] + + combo_info: list[dict] = [] # 每个组合的生成信息(供 F2/BC 亲本 BV、系谱/区试映射) + promo_queue: list[dict] = [] # 正向晋级队列(按组合内轮次顺序,最终走服务层) + elim_accum: list[dict] = [] # 淘汰记录(全部晋级完成后统一批量回写) + + # ---------- F1 ---------- + f1_combos = [(y, seq) for y in f1_years for seq in range(1, f1_per_year + 1)] + print(f"[sim] 3/8 生成 F1 杂交组合 {len(f1_years)} 年度 × {f1_per_year} 组合/年 = {len(f1_combos)} 个 → 单株 ...") + combo_rows = [] + for y, seq in f1_combos: + dam_i = (y * 7 + seq * 3) % nf + sire_i = (dam_i + 1 + seq) % nf + while sire_i == dam_i: + sire_i = (sire_i + 1) % nf + combo_rows.append({ + "combination_code": f"桃{y}-{seq:03d}", "cross_year": y, + "bre_target_id": target_ids[(y + seq) % len(target_ids)], + "female_parent_id": founder_ids[dam_i], "male_parent_id": founder_ids[sire_i], + "cross_method": "人工授粉", "cross_type": "杂交", "design_type": "full_diallel", + "stage": "seedling", "cross_date": f"{y}-04-{10 + seq % 10:02d}", + "seed_count": pyrand.randint(300, 800), + "reason": "目标亲本组配", "remark": "F1 杂种实生苗圃", + }) + res = await session.execute( + insert(CrossCombinationModel).values([{**_base(), **r} for r in combo_rows]).returning(CrossCombinationModel.id)) + f1_combo_ids = list(res.scalars()) + counts["combos"] += len(f1_combo_ids) + + combo_defs = [] + for (y, seq), cid in zip(f1_combos, f1_combo_ids): + code = f"桃{y}-{seq:03d}" + dam_i = (y * 7 + seq * 3) % nf + sire_i = (dam_i + 1 + seq) % nf + while sire_i == dam_i: + sire_i = (sire_i + 1) % nf + combo_defs.append({"code": code, "id": cid, "generation": "F1", + "dam_i": dam_i, "sire_i": sire_i, + "dam_id": founder_ids[dam_i], "sire_id": founder_ids[sire_i], + "planted": y + 2, "cross_year": y, "n": f1_trees}) + + for c in combo_defs: + await _gen_combo(session, counts, combo_info, promo_queue, elim_accum, c, + founder_bv, params, sigmaA, sd, w, codes, rng, pyrand, + personnel_ids, plots, rule_ids, rootstock_ids, end_year, trait_id) + await session.commit() + print(f"[sim] 累计: 组合={counts['combos']} 单株={counts['trees']} 评价={counts['evals']} " + f"观测={counts['obs']} 农事={counts['field_ops']} 淘汰={len(elim_accum)}") + + # ---------- F1 选择晋级(真实服务层 create_batch) ---------- + print("[sim] 4/8 运行 F1 选择晋级(服务层批量,sp/ap/line/regional/released)...") + ok, fail = await run_promotions(session, auth, promo_queue, args.batch_size) + counts["selections"] += ok + promo_queue.clear() + print(f"[sim] F1 晋级成功 {ok} / 失败 {fail}") + + # 收集 line 及以上种质(id / BV),供 F2/BC 亲本与区试 + line_germ = await collect_line_germplasms(session, combo_info) + print(f"[sim] line 及以上种质 {len(line_germ)} 个") + + # ---------- F2 / BC ---------- + print("[sim] 5/8 生成 F2 自交 / BC1 回交组合与单株 ...") + d2c = await _gen_derived_combos( + session, counts, combo_info, promo_queue, elim_accum, line_germ, founder_bv, founder_ids, + target_ids, f2_combos, f2_trees, bc_combos, bc_trees, + params, sigmaA, sd, w, codes, rng, pyrand, personnel_ids, plots, + rule_ids, rootstock_ids, end_year, trait_id) + await session.commit() + print(f"[sim] 派生组合 {d2c} 个,累计单株 {counts['trees']} 观测 {counts['obs']}") + + print("[sim] 6/8 运行 F2/BC 选择晋级 ...") + ok2, fail2 = await run_promotions(session, auth, promo_queue, args.batch_size) + counts["selections"] += ok2 + promo_queue.clear() + print(f"[sim] F2/BC 晋级成功 {ok2} / 失败 {fail2}") + + # 晋级全部完成后:批量回写淘汰记录 + 淘汰树状态 + if elim_accum: + await insert_chunked(session, SelectionResultModel, [{**_base(), **r} for r in elim_accum]) + counts["selections"] += len(elim_accum) + await apply_elim_status(session) + await session.commit() + print(f"[sim] 淘汰记录 {len(elim_accum)} 条已回写,淘汰树状态已置位") + + # ---------- 区域试验 ---------- + print("[sim] 7/8 生成区域试验(trial / study / entry)...") + await gen_trials(session, counts, combo_info, plots, target_ids, plan_ids) + await session.commit() + + # ---------- 业务流空表补齐 ---------- + print("[sim] 8/8 补齐业务流空表(授粉/花粉/种子批/种子处理/育苗/克隆扩繁/环境气象/通用观测/试验处理)...") + await gen_support_data(session, counts, combo_info, personnel_ids, plots, + rootstock_ids, trait_id, params, codes, end_year, pyrand, rng) + await session.commit() + + # ---------- 自检 ---------- + elapsed = (datetime.now(UTC) - t0).total_seconds() + print(f"[sim] 生成完成,耗时 {elapsed:.1f}s。") + await report_counts(session, counts) + if args.verify: + await verify(session, counts, dry) + + await engine.dispose() + + +# --------------------------------------------------------------------------- +# 单个组合生成(组合 → 定植 → 单株 → 选择排程 → 淘汰记录/评价/观测/农事 + 晋级入队) +# --------------------------------------------------------------------------- +async def _gen_combo(session, counts, combo_info, promo_queue, elim_accum, c, + founder_bv, params, sigmaA, sd, w, codes, rng, pyrand, + personnel_ids, plots, rule_ids, rootstock_ids, end_year, trait_id) -> None: + code = c["code"] + cid = c["id"] + n = c["n"] + planted = c["planted"] + generation = c["generation"] + + # 父本 BV + if generation == "F1": + dam_bv, sire_bv = founder_bv[c["dam_i"]], founder_bv[c["sire_i"]] + else: + dam_bv, sire_bv = c["dam_bv"], c["sire_bv"] + + # 单株 BV 与指数 + bv = np.zeros((n, len(codes))) + for i in range(n): + bv[i] = sample_bv_parents(rng, dam_bv, sire_bv, sigmaA) + idx = idx_score(bv, sd, w) + + plot = plots[c["cross_year"] % len(plots)] + root = rootstock_ids[(c["cross_year"] * 3) % len(rootstock_ids)] + persona = personnel_ids[c["cross_year"] % len(personnel_ids)] + + # 定植批次 + plant_row = {**_base(), + "combination_id": cid, "plot_id": plot["id"], "planting_date": f"{planted}-03-15", + "tree_count": n, "row_no": 1, "col_no": 1, "bre_personnel_id": persona, + "rootstock_id": root, "block_no": 1, "remark": f"{code} 定植"} + pres = await session.execute(insert(PlantingModel).values([plant_row]).returning(PlantingModel.id)) + planting_id = list(pres.scalars())[0] + counts["plantings"] += 1 + + # 单株 + tree_rows = [{ + "combination_id": cid, "dam_id": c["dam_id"], "sire_id": c["sire_id"], + "tree_no": f"{code}-{i + 1:03d}", "plot_id": plot["id"], "planting_id": planting_id, + "row_no": i // 20 + 1, "col_no": i % 20 + 1, "rootstock_id": root, + "planted_date": f"{planted}-03-15", "bre_personnel_id": persona, + "status": "alive", "stage": "seedling", "generation": generation, + "remark": f"{code} 实生苗", + } for i in range(n)] + res = await session.execute(insert(TreeModel).values([{**_base(), **r} for r in tree_rows]).returning(TreeModel.id)) + tree_ids = list(res.scalars()) + counts["trees"] += n + + info = {"code": code, "combo_id": cid, "generation": generation, "planted": planted, + "cross_year": c["cross_year"], "n": n, + "tree_ids": tree_ids, "bv": bv, "rounds": {}} + combo_info.append(info) + + # 童期观测(定植年 ~ sp 前一年 / 窗口末) + sp_year = planted + 3 + await gen_juvenile_obs(session, counts, cid, tree_ids, bv, trait_id, params, sd, w, + planted, min(sp_year - 1, end_year), rng, code) + + # 选择排程 + alive_pool = np.arange(n) + current_sel = None + for ri, (rname, is_selected, frac, from_stage, to_stage, off) in enumerate(ROUNDS): + ryear = planted + off + pool = alive_pool if rname == "sp" else current_sel + if pool is None or ryear > end_year: + info["rounds"][rname] = {"year": ryear, "active": False} + continue + sel, elim = select_pool(pool, idx, frac, rng) + # 淘汰记录(暂存,晋级全部完成后统一回写) + for ei in elim: + elim_accum.append({ + "combination_id": cid, "tree_id": tree_ids[ei], "tree_no": f"{code}-{ei + 1:03d}", + "selection_year": ryear, "is_selected": "eliminated", + "from_stage": from_stage, "to_stage": from_stage, + "reason": f"{rname} 轮未达入选标准", "remark": "模拟数据"}) + info["rounds"][rname] = {"year": ryear, "active": True, "sel": sel, "elim": elim} + if rname == "sp": + info["sp_order"] = [int(i) for i in sel] # 晋级 clone 序由此决定 + # 入队正向晋级(选中树) + if len(sel): + promo_queue.append({ + "combination_id": cid, "selection_year": ryear, "is_selected": is_selected, + "from_stage": from_stage, "to_stage": to_stage, + "rule_id": rule_ids[min(ri, len(rule_ids) - 1)] if rule_ids else None, + "approved_by": personnel_ids[0], + "reason": f"{rname} 轮入选(指数前 {len(sel)} 名)", "remark": "模拟数据", + "tree_ids": [tree_ids[int(i)] for i in sel], + }) + current_sel = sel + + # 成株评价与观测:各轮选中树,自 sp_year 至淘汰年前一年/窗口末 + await gen_eval_obs(session, counts, info, cid, tree_ids, bv, trait_id, params, sd, w, + sp_year, end_year, rng, personnel_ids, plot, code, pyrand) + + # 小区级农事(每年) + await gen_plot_ops(session, counts, plot["id"], planted, end_year, personnel_ids, pyrand, code) + + +def _gen_desc_obs_rows(base, tree_id, cid, y, trait_id, code, bvec, idx_of, + pyrand, rng, year_eff=None, reuse_vals=None, + evaluation_id=None, stage="evaluation") -> list[dict]: + """为单棵树单年份生成 43 条描述性状观测(仅入选候选树 / ramet 分层使用)。 + + categorical → value_text(scale_json 选项抽样);date → value_date(花期 3-4 月); + numeric core → value_numeric(遗传代理 + 年度 + 误差);纵/横径复用 stats 同源观测值。 + """ + rows = [] + year_eff = year_eff or {} + reuse_vals = reuse_vals or {} + for dcode, _name, _cat, dtype, options, mu, sd, vmin, vmax, reuse in DESCRIPTIVE_TRAITS: + # avg_fruit_weight/ssc 与 stats 撞 code 且同为 numeric:其观测由既有数值循环写入 + # (同一 (evaluation_id, trait_id) 唯一约束下不可双写),描述侧不再造值。 + if dcode in ("avg_fruit_weight", "ssc"): + continue + row = {**_base(), + "tree_id": tree_id, "combination_id": cid, + "trait_id": trait_id[dcode], "evaluate_year": y, + "stage": stage, "remark": f"{code} 描述观测"} + if evaluation_id is not None: + row["evaluation_id"] = evaluation_id + if dtype == "categorical": + row["value_text"] = str(pyrand.choice(options)) + elif dtype == "date": + row["value_date"] = f"{y}-{pyrand.randint(3, 4):02d}-{pyrand.randint(1, 28):02d}" + else: + if reuse is not None: + v = reuse_vals.get(reuse) + if v is None: + continue + row["value_numeric"] = float(v) + else: + proxy, scale = DESCRIPTIVE_BV_PROXY.get(dcode, (None, 1.0)) + g = scale * bvec[idx_of[proxy]] if proxy is not None and proxy in idx_of else 0.0 + row["value_numeric"] = clip_round(mu + g + year_eff.get(dcode, 0.0) + + rng.normal(0, sd * 0.6), vmin, vmax) + rows.append(row) + return rows + + +async def gen_juvenile_obs(session, counts, cid, tree_ids, bv, trait_id, params, sd, w, + start_year, end_year, rng, code) -> None: + if end_year < start_year: + return + juv = [c for c in params if params[c]["stage"] == "juvenile"] + rows = [] + for y in range(start_year, end_year + 1): + year_eff = {c: rng.normal(0, params[c]["sd"] * 0.10) for c in juv} + for i in range(len(tree_ids)): + for c in juv: + p = params[c] + val = clip_round(p["mu"] + bv[i, p["idx"]] + year_eff[c] + rng.normal(0, p["sd"] * 0.6), p["vmin"], p["vmax"]) + rows.append({**_base(), + "tree_id": tree_ids[i], "combination_id": cid, + "trait_id": trait_id[c], "evaluate_year": y, + "value_numeric": val, "stage": "juvenile", + "remark": f"{code} 童期观测"}) + if rows: + counts["obs"] += await insert_chunked(session, TraitObservationModel, rows) + + +async def gen_eval_obs(session, counts, info, cid, tree_ids, bv, trait_id, params, sd, w, + sp_year, end_year, rng, personnel_ids, plot, code, pyrand) -> None: + eval_num = [c for c in params if params[c]["stage"] == "evaluation"] + text_specs = {t[0]: t for t in TEXT_TRAITS} + idx_of = {c: params[c]["idx"] for c in params} + eval_tree_idx = None + for rname in ("sp", "ap", "line", "regional", "released"): + r = info["rounds"].get(rname) + if r and r.get("active"): + eval_tree_idx = r["sel"] + break + if eval_tree_idx is None: + return + # 各树淘汰年前一年为止(无淘汰则评估至窗口末) + elim_year = {} + for rname in ("ap", "line", "regional", "released"): + r = info["rounds"].get(rname) + if r and r.get("active") and r.get("elim") is not None: + for ei in r["elim"]: + elim_year[ei] = r["year"] - 1 + for y in range(sp_year, end_year + 1): + active = [i for i in eval_tree_idx if elim_year.get(i, end_year + 1) >= y] + if not active: + continue + year_eff = {c: rng.normal(0, params[c]["sd"] * 0.10) for c in eval_num} + # 单株评价(先插,取 id 供观测外键) + e_rows, eval_trees = [], [] + for i in active: + z = idx_score(bv[i:i + 1], sd, w)[0] + score = float(np.clip(55 + 25 * (z / 1.5) + rng.normal(0, 3), 0, 100)) + m = pyrand.randint(6, 8) + d = pyrand.randint(10, 28) + e_rows.append({**_base(), + "combination_id": cid, "tree_id": tree_ids[i], "evaluate_year": y, + "evaluate_date": f"{y}-{m:02d}-{d:02d}", + "bre_personnel_id": personnel_ids[y % len(personnel_ids)], + "overall_score": round(score, 1), + "crop_load": round(pyrand.uniform(0.3, 1.0), 2), + "remark": f"{code} 成株评价"}) + eval_trees.append((i, m, d)) + res = await session.execute(insert(TreeEvaluationModel).values(e_rows).returning(TreeEvaluationModel.id)) + eval_ids = list(res.scalars()) + counts["evals"] += len(e_rows) + eval_map = {i: eid for (i, _m, _d), eid in zip(eval_trees, eval_ids)} + date_map = {i: (m, d) for i, m, d in eval_trees} + # 观测(数值 + 文本)与采收(同源表型,保证一致) + o_rows, h_rows = [], [] + for i in active: + eid = eval_map[i] + vals = {} + for c in eval_num: + p = params[c] + val = clip_round(p["mu"] + bv[i, p["idx"]] + year_eff[c] + rng.normal(0, p["sd"] * 0.6), p["vmin"], p["vmax"]) + vals[c] = val + o_rows.append({**_base(), + "evaluation_id": eid, "tree_id": tree_ids[i], "combination_id": cid, + "trait_id": trait_id[c], "evaluate_year": y, + "value_numeric": val, "stage": "evaluation", + "crop_load": round(pyrand.uniform(0.3, 1.0), 2), + "remark": f"{code} 成株观测"}) + for tcode, _name, _cat, _st, opts, probs in text_specs.values(): + if tcode in ("flesh_color", "fruit_shape"): + continue # 已由描述性状 categorical 覆盖(同一 (evaluation_id, trait_id) 不双写) + o_rows.append({**_base(), + "evaluation_id": eid, "tree_id": tree_ids[i], "combination_id": cid, + "trait_id": trait_id[tcode], "evaluate_year": y, + "value_text": str(pyrand.choices(opts, weights=probs, k=1)[0]), + "stage": "evaluation", "remark": f"{code} 成株观测"}) + # 描述性状:候选树全量 43 条,挂同一评价头(纵/横径复用本次 stats 观测) + o_rows.extend(_gen_desc_obs_rows(_base(), tree_ids[i], cid, y, trait_id, code, + bv[i], idx_of, pyrand, rng, + year_eff=year_eff, reuse_vals=vals, + evaluation_id=eid, stage="evaluation")) + m, d = date_map[i] + h_rows.append({**_base(), + "tree_id": tree_ids[i], "plot_id": plot["id"], "op_type": "harvest", + "op_date": date(y, m, d), "op_detail": "采收并记录单株产量", + "yield_kg": vals.get("single_tree_yield"), + "fruit_count": int(vals.get("fruit_number") or 0), + "avg_fruit_weight": vals.get("avg_fruit_weight"), + "marketable_rate": vals.get("marketable_rate"), + "operator_id": personnel_ids[y % len(personnel_ids)], "status": 1, + "remark": f"{code} 采收"}) + if o_rows: + counts["obs"] += await insert_chunked(session, TraitObservationModel, o_rows) + if h_rows: + counts["field_ops"] += await insert_chunked(session, FieldOperationModel, h_rows) + + +async def gen_plot_ops(session, counts, plot_id, start_year, end_year, personnel_ids, pyrand, code) -> None: + ops = [("fertilize", 3), ("spray", 5), ("prune", 2), ("irrigate", 6)] + rows = [] + for y in range(start_year, end_year + 1): + for op, mo in ops: + rows.append({**_base(), "plot_id": plot_id, "op_type": op, + "op_date": date(y, mo, pyrand.randint(1, 28)), + "op_detail": f"{code} 地块 {op}", "operator_id": personnel_ids[y % len(personnel_ids)], + "status": 1, "remark": "小区例行农事"}) + if rows: + counts["field_ops"] += await insert_chunked(session, FieldOperationModel, rows) + + +# --------------------------------------------------------------------------- +# 批量晋级(真实服务层 create_batch,事务化 + 失败单株回退重试) +# --------------------------------------------------------------------------- +async def run_promotions(session, auth, queue, batch_size) -> tuple[int, int]: + ok = fail = 0 + with bre_audit_suppress(): + for item in queue: + base = {k: item[k] for k in + ("selection_year", "is_selected", "rule_id", "approved_by", + "from_stage", "to_stage", "reason", "remark")} + for chunk in chunked(item["tree_ids"], batch_size): + data = SelectionResultBatchCreateSchema(combination_id=item["combination_id"], tree_id=chunk, **base) + try: + svc = SelectionResultService(auth, session) + res = await svc.create_batch(data) + await session.commit() + ok += res.success_count + fail += res.fail_count + for fd in res.fail_details: + print(f"[sim] !! 晋级失败: {fd}") + except Exception as e: # 整批失败:逐棵回退重试 + await session.rollback() + for tid in chunk: + single = SelectionResultBatchCreateSchema(combination_id=item["combination_id"], tree_id=[tid], **base) + try: + svc = SelectionResultService(auth, session) + r = await svc.create_batch(single) + await session.commit() + ok += r.success_count + fail += r.fail_count + except Exception as e2: + await session.rollback() + fail += 1 + print(f"[sim] !! 单棵晋级失败 tree_id={tid}: {e2}") + return ok, fail + + +# --------------------------------------------------------------------------- +# 派生组合(F2 / BC) +# --------------------------------------------------------------------------- +async def collect_line_germplasms(session, combo_info) -> list[dict]: + """返回 line 及以上选中树的种质信息 [{code, seq, bv, germplasm_id, combo_id, year, clone_code}]。""" + names = [] + for info in combo_info: + r = info["rounds"].get("line") + if not r or not r.get("active") or "sel" not in r: + continue + sp_order = info.get("sp_order", []) + for ei in r["sel"]: + ei = int(ei) + if ei not in sp_order: + continue + seq = sp_order.index(ei) + 1 + names.append({"code": info["code"], "seq": seq, "bv": info["bv"][ei].copy(), + "combo_id": info["combo_id"], "year": r["year"], "clone_code": f"{info['code']}-{seq:04d}"}) + if not names: + return [] + rows = await session.execute( + text("SELECT id, cultivar_name FROM bre_germplasm WHERE cultivar_name = ANY(:codes)"), + {"codes": [n["clone_code"] for n in names]}) + gid_by_name = {name: gid for gid, name in rows.all()} + for n in names: + n["germplasm_id"] = gid_by_name.get(n["clone_code"]) + return [n for n in names if n["germplasm_id"]] + + +async def _gen_derived_combos(session, counts, combo_info, promo_queue, elim_accum, line_germ, + founder_bv, founder_ids, target_ids, f2_combos, f2_trees, + bc_combos, bc_trees, + params, sigmaA, sd, w, codes, rng, pyrand, personnel_ids, plots, + rule_ids, rootstock_ids, end_year, trait_id) -> int: + if not line_germ: + return 0 + candidates = [g for g in line_germ if g["year"] <= end_year - 3] + if not candidates: + return 0 + nfg = len(founder_ids) + ncan = len(candidates) + defs = [] + for k in range(f2_combos): + a = candidates[k % ncan] + b_idx = (k * 7 + 3) % ncan + b = candidates[b_idx] + # 候选池中同一组合的多个 line 种质(同 combo_id)连续成块, + # 若固定 index 取到同块元素会无限循环,故循环体内 index 递增并加防死循环上限 + guard = 0 + while b["combo_id"] == a["combo_id"] and guard < ncan: + b_idx = (b_idx + 1) % ncan + b = candidates[b_idx] + guard += 1 + if b["combo_id"] == a["combo_id"]: + continue + cy = max(a["year"], b["year"]) + 1 + if cy > end_year - 1: + continue + defs.append({"generation": "F2", "dam": a, "sire": b, "cross_year": cy, "n": f2_trees, "kind": "F2"}) + for k in range(bc_combos): + a = candidates[k % len(candidates)] + fi = (k * 5 + 1) % nfg + cy = a["year"] + 1 + if cy > end_year - 1: + continue + defs.append({"generation": "BC1", "dam": a, + "sire": {"bv": founder_bv[fi], "germplasm_id": founder_ids[fi]}, + "cross_year": cy, "n": bc_trees, "kind": "BC"}) + if not defs: + return 0 + combo_rows = [] + for k, d in enumerate(defs): + code = f"桃{d['cross_year']}{d['kind']}-{k + 1:03d}" + combo_rows.append({ + "combination_code": code, "cross_year": d["cross_year"], + "bre_target_id": target_ids[k % len(target_ids)], + "female_parent_id": d["dam"].get("germplasm_id"), + "male_parent_id": d["sire"].get("germplasm_id"), + "cross_method": "人工授粉", "cross_type": "杂交" if d["kind"] == "F2" else "回交", + "design_type": "partial_diallel", "stage": "seedling", + "cross_date": f"{d['cross_year']}-04-{10 + k % 10:02d}", "seed_count": pyrand.randint(200, 500), + "reason": "世代推进", "remark": f"{d['kind']} 分离世代", + }) + res = await session.execute( + insert(CrossCombinationModel).values([{**_base(), **r} for r in combo_rows]).returning(CrossCombinationModel.id)) + ids = list(res.scalars()) + counts["combos"] += len(ids) + nd = len(defs) + for k, (d, cid) in enumerate(zip(defs, ids)): + if k % 20 == 0 or k == nd - 1: + print(f"[sim] F2/BC 组合 {k + 1}/{nd}({d['kind']} {d['cross_year']})...") + c = {"code": f"桃{d['cross_year']}{d['kind']}-{k + 1:03d}", "id": cid, + "generation": d["generation"], "planted": d["cross_year"] + 1, + "dam_bv": d["dam"]["bv"], "sire_bv": d["sire"]["bv"], + "dam_id": d["dam"].get("germplasm_id"), "sire_id": d["sire"].get("germplasm_id"), + "cross_year": d["cross_year"], "n": d["n"]} + await _gen_combo(session, counts, combo_info, promo_queue, elim_accum, c, + founder_bv, params, sigmaA, sd, w, codes, rng, pyrand, + personnel_ids, plots, rule_ids, rootstock_ids, end_year, trait_id) + return len(ids) + + +# --------------------------------------------------------------------------- +# 区域试验 +# --------------------------------------------------------------------------- +async def gen_trials(session, counts, combo_info, plots, target_ids, plan_ids) -> None: + regional = [] # (combo_id, year, clone_code) + for info in combo_info: + r = info["rounds"].get("regional") + if not r or not r.get("active") or "sel" not in r: + continue + sp_order = info.get("sp_order", []) + for ei in r["sel"]: + ei = int(ei) + if ei not in sp_order: + continue + seq = sp_order.index(ei) + 1 + regional.append((info["combo_id"], r["year"], f"{info['code']}-{seq:04d}")) + if not regional: + return + rows = await session.execute( + text("SELECT id, cultivar_name FROM bre_germplasm WHERE cultivar_name = ANY(:codes)"), + {"codes": [name for _, _, name in regional]}) + gid_by_name = {name: gid for gid, name in rows.all()} + by_year: dict[int, list] = {} + for cid, y, name in regional: + gid = gid_by_name.get(name) + if gid: + by_year.setdefault(y, []).append((cid, gid)) + for y, items in sorted(by_year.items()): + trial = {**_base(), "trial_name": f"桃区试-{y}", "plan_id": plan_ids[0] if plan_ids else None, + "target_id": target_ids[y % len(target_ids)] if target_ids else None, + "trial_type": "regional", "design_type": "randomized_block", + "start_year": y, "end_year": y + 2, "objective": f"{y} 年区域试验", "remark": "模拟数据"} + tres = await session.execute(insert(TrialModel).values([trial]).returning(TrialModel.id)) + trial_id = list(tres.scalars())[0] + counts["trials"] += 1 + for site in {p["site_id"] for p in plots}: + study = {**_base(), "trial_id": trial_id, "study_name": f"区试-{y}-site{site}", "site_id": site, + "year": y, "block_count": 3, "season": "夏", "design_type": "randomized_block", + "remark": "模拟数据"} + sres = await session.execute(insert(TrialStudyModel).values([study]).returning(TrialStudyModel.id)) + study_id = list(sres.scalars())[0] + counts["studies"] += 1 + entry_rows = [{**_base(), "trial_study_id": study_id, "entry_number": eno, + "germplasm_id": gid, "combination_id": cid, "block_no": eno % 3 + 1, + "remark": "模拟数据"} + for eno, (cid, gid) in enumerate(items, start=1)] + await insert_chunked(session, TrialStudyEntryModel, entry_rows) + counts["entries"] += len(entry_rows) + + +# --------------------------------------------------------------------------- +# 业务流空表补齐(授粉→花粉→种子批→种子处理→育苗、克隆扩繁、环境气象、通用观测、试验处理) +# --------------------------------------------------------------------------- +async def gen_support_data(session, counts, combo_info, personnel_ids, plots, rootstock_ids, + trait_id, params, codes, end_year, pyrand, rng) -> None: + """为流程侧空表生成模拟数据(纯批量 insert,不走服务层)。 + + 授粉链按组合逐条生成(F1/F2/BC 全含);扩繁只对晋级到 line 的树(其 clone 已由 + 选择晋级服务层创建);气象按 site×year 唯一;通用观测给 line 树补 EAV 记录。 + """ + # 组合基础信息(seed_count / 父本种质) + combo_meta = {cid: {"seed_count": sc, "male_parent_id": mpid} + for cid, sc, mpid in (await session.execute( + text("SELECT id, seed_count, male_parent_id FROM bre_cross_combination"))).all()} + name_by_gid = {gid: name for gid, name in (await session.execute( + select(BreedingGermplasmModel.id, BreedingGermplasmModel.cultivar_name))).all()} + + # 1) 每组合一条授粉链:pollen → pollination → seed_lot → seed_treatment → seedling + pollen_rows, polli_rows, seed_rows, treat_rows, seedl_rows = [], [], [], [], [] + for info in combo_info: + code, cid = info["code"], info["combo_id"] + cy = info["cross_year"] + meta = combo_meta.get(cid, {"seed_count": pyrand.randint(300, 800), "male_parent_id": None}) + seed_count = meta["seed_count"] or pyrand.randint(300, 800) + mpid = meta.get("male_parent_id") + plot = plots[cy % len(plots)] + persona = personnel_ids[cy % len(personnel_ids)] + # 花粉批次(父本混合花粉,母/父本为种质无树) + pollen_rows.append({**_base(), + "lot_code": f"PL-{code}", "male_tree_id": None, "source_type": "mixed", + "source_desc": f"父本混合花粉({name_by_gid.get(mpid, '亲本')})", + "collect_date": date(cy, 3, 15), "collect_method": "采花取粉", + "quantity": f"{pyrand.randint(20, 80)}g", "storage_method": "minus80", + "viability_method": "ttc", "viability_pct": round(pyrand.uniform(70, 92), 1), + "viability_test_date": date(cy, 3, 20), "expiry_date": date(cy, 5, 15), + "remark": "模拟数据"}) + # 授粉 + polli_rows.append({**_base(), + "combination_id": cid, "plot_id": plot["id"], + "pollination_date": f"{cy}-04-{10 + cy % 10:02d}", + "flower_count": pyrand.randint(800, 3000), "effective_count": pyrand.randint(200, 800), + "bre_personnel_id": persona, "pollinator": None, + "emasculation_date": date(cy, 4, 5), "bagging_date": date(cy, 4, 25), + "female_tree_id": None, "male_tree_id": None, + "pollination_method": "人工点授", "pollen_lot_id": None, "remark": "模拟数据"}) + # 育苗出苗数(衔接定植株数 n),种子批 used_count 同步 + emerg = max(int(info["n"] * pyrand.uniform(1.2, 2.0)), 1) + strong = int(emerg * pyrand.uniform(0.85, 1.0)) + seed_rows.append({**_base(), + "combination_id": cid, "pollination_id": None, "lot_code": f"SL-{code}", + "harvest_year": cy, "seed_count": seed_count, "used_count": emerg, + "germination_rate": round(pyrand.uniform(55, 85), 2), + "storage_type": "种子库", "storage_location": f"库A-{cy % 8 + 1}架", + "test_date": date(cy, 8, 15), "remark": "模拟数据"}) + # 种子处理(冬季低温层积) + treat_rows.append({**_base(), + "combination_id": cid, "seed_lot_id": None, "treatment_method": "低温层积催芽", + "treatment_start": f"{cy}-10-01", "treatment_end": f"{cy + 1}-02-15", + "germination_rate": round(pyrand.uniform(55, 85), 1), + "bre_personnel_id": persona, "remark": "模拟数据"}) + # 育苗(stage 3 壮苗,衔接定植) + sowing = f"{cy + 1}-03-05" + seedl_rows.append({**_base(), + "combination_id": cid, "treatment_id": None, "seed_lot_id": None, + "batch_no": f"YP-{cid}-{sowing.replace('-', '')}-1", + "sowing_date": sowing, "tray_no": f"T-{cy % 40 + 1:03d}", "nursery": "育苗圃", + "seedling_count": emerg, "emergence_date": f"{cy + 1}-04-01", + "strong_seedling_date": f"{cy + 1}-05-10", "strong_seedling_count": strong, + "stage": "3", "bre_personnel_id": persona, "remark": "模拟数据"}) + + if pollen_rows: + pollen_ids = list((await session.execute( + insert(PollenModel).values(pollen_rows).returning(PollenModel.id))).scalars()) + for r, pid in zip(polli_rows, pollen_ids): + r["pollen_lot_id"] = pid + polli_ids = list((await session.execute( + insert(PollinationModel).values(polli_rows).returning(PollinationModel.id))).scalars()) + for r, pid in zip(seed_rows, polli_ids): + r["pollination_id"] = pid + seed_ids = list((await session.execute( + insert(SeedLotModel).values(seed_rows).returning(SeedLotModel.id))).scalars()) + for r, sid in zip(treat_rows, seed_ids): + r["seed_lot_id"] = sid + treat_ids = list((await session.execute( + insert(SeedTreatmentModel).values(treat_rows).returning(SeedTreatmentModel.id))).scalars()) + for r, tid, sid in zip(seedl_rows, treat_ids, seed_ids): + r["treatment_id"], r["seed_lot_id"] = tid, sid + await insert_chunked(session, SeedlingModel, seedl_rows) + counts.update(pollen=len(pollen_ids), pollination=len(polli_ids), + seed_lots=len(seed_ids), seed_treatments=len(treat_ids), + seedlings=len(seedl_rows)) + await session.commit() + + # 2) 克隆扩繁闭环:为每个入选克隆补 K 棵共享 clone_id 的无性系苗(ramet)树 + 扩繁批次 + # + 童期/成株观测。ABLUP G1 门禁按 c{clone}(回退 f{组合})分组且要求每组 ≥ min_clone_n, + # 服务层晋级为每入选树各建唯一克隆(1:1),缺 ramet 即每组 n=1 无法建模,故按真实 + # "入选株→克隆→扩繁→新树"闭环补齐;ramet 共享克隆 BV,观测仅加年度与误差噪声。 + clone_id_by_code = {cc: cid for cid, cc in (await session.execute( + select(CloneModel.id, CloneModel.clone_code))).all()} + orig_by_clone: dict[int, tuple] = {} + for tid, cl, cid_, gen, dam, sire, plot_id, root, persona, stg in (await session.execute( + text("SELECT id, clone_id, combination_id, generation, dam_id, sire_id, plot_id, " + "rootstock_id, bre_personnel_id, stage FROM bre_tree " + "WHERE clone_id IS NOT NULL AND is_deleted = false"))).all(): + orig_by_clone[cl] = (tid, cid_, gen, dam, sire, plot_id, root, persona, stg) + eval_num = [c for c in params if params[c]["stage"] == "evaluation"] + juv_num = [c for c in params if params[c]["stage"] == "juvenile"] + idx_of = {c: params[c]["idx"] for c in params} + prop_rows = [] + for info in combo_info: + sp_order = info.get("sp_order", []) + if not sp_order: + continue + code, cy = info["code"], info["cross_year"] + sp_year = info["planted"] + 3 + ev_years = [y for y in range(sp_year, min(sp_year + RAMET_EVAL_YEARS, end_year + 1))] + if not ev_years: + continue + bv = info["bv"] + combo_id = info["combo_id"] + gen = info["generation"] + plot = plots[cy % len(plots)] + local_ramets, local_meta = [], [] + for seq, ei in enumerate(sp_order, 1): + clone_code = f"{code}-{seq:04d}" + cid_ = clone_id_by_code.get(clone_code) + orig = orig_by_clone.get(cid_) if cid_ else None + if not orig: + continue + _tid, _cid, _gen, dam, sire, plot_id, root, persona, stg = orig + rs_id = root or rootstock_ids[(cy + seq) % len(rootstock_ids)] + op_id = persona or personnel_ids[(cy + seq) % len(personnel_ids)] + grafted = pyrand.randint(40, 90) + prop_rows.append({**_base(), + "batch_code": f"PR-{clone_code}", "scion_source_type": "tree", + "scion_source_id": _tid, "produced_clone_id": cid_, + "rootstock_id": rs_id, "method": "嫁接", "graft_date": date(sp_year, 3, 5), + "nursery_site_id": plot["site_id"], "operator_id": op_id, + "scion_count": grafted, "grafted_count": grafted, + "survival_count": int(grafted * pyrand.uniform(0.75, 0.95)), + "destination": "定植", "remark": "模拟数据"}) + for k in range(K_RAMETS): + local_ramets.append({**_base(), + "combination_id": combo_id, "dam_id": dam, "sire_id": sire, + "clone_id": cid_, "plot_id": plot_id or plot["id"], + "tree_no": f"{code}-{seq:04d}R{k + 1:02d}", + "row_no": k + 1, "col_no": seq % 20 + 1, "block_no": seq % 3 + 1, + "rootstock_id": rs_id, "planted_date": f"{sp_year}-03-10", + "bre_personnel_id": op_id, "status": "alive", "stage": stg, + "generation": gen, "remark": f"{clone_code} 克隆无性系苗"}) + local_meta.append((bv[ei], ev_years, code, combo_id)) + if prop_rows: + await insert_chunked(session, PropagationModel, prop_rows) + counts["propagations"] += len(prop_rows) + prop_rows = [] + if local_ramets: + rid_rows = list((await session.execute( + insert(TreeModel).values(local_ramets).returning(TreeModel.id))).scalars()) + counts["trees"] += len(rid_rows) + obs = [] + for tid, (bvec, years, code_, cid_) in zip(rid_rows, local_meta): + for y in years: + year_eff = {c: rng.normal(0, params[c]["sd"] * 0.10) + for c in eval_num + juv_num} + reuse_vals = {} + for c in eval_num + juv_num: + p = params[c] + val = clip_round(p["mu"] + bvec[p["idx"]] + year_eff[c] + + rng.normal(0, p["sd"] * 0.6), p["vmin"], p["vmax"]) + if c in ("fruit_length", "fruit_diameter"): + reuse_vals[c] = val + obs.append({**_base(), + "tree_id": tid, "combination_id": cid_, "trait_id": trait_id[c], + "evaluate_year": y, "value_numeric": val, + "stage": params[c]["stage"], "remark": f"{code_} 克隆观测"}) + # 描述性状:克隆苗全量 43 条,无评价头(tree_id 直挂),纵/横径复用本次观测 + obs.extend(_gen_desc_obs_rows(_base(), tid, cid_, y, trait_id, code_, + bvec, idx_of, pyrand, rng, + year_eff=year_eff, reuse_vals=reuse_vals, + stage="evaluation")) + counts["obs"] += await insert_chunked(session, TraitObservationModel, obs, 20000) + + # 3) 环境气象:site × year(2007..end_year,site×year 唯一) + env_rows = [] + for sid in sorted({p["site_id"] for p in plots}): + for y in range(2007, end_year + 1): + env_rows.append({**_base(), + "site_id": sid, "year": y, + "chilling_hours": round(750 + sid * 40 + pyrand.uniform(-120, 120), 1), + "growing_degree_days": round(3200 + pyrand.uniform(-300, 300), 1), + "rainfall_mm": round(620 + sid * 80 + pyrand.uniform(-120, 120), 1), + "temp_avg": round(14.5 + sid * 0.8 + pyrand.uniform(-1.5, 1.5), 1), + "soil_moisture": round(0.55 + pyrand.uniform(-0.1, 0.1), 2), + "source": "气象站插值", "remark": "模拟数据"}) + if env_rows: + await insert_chunked(session, EnvironmentConditionModel, env_rows) + counts["env_conds"] += len(env_rows) + + # 4) 通用观测(EAV):line 晋级树补 ssc(数值) / 果肉颜色(文本) / 果实形状(文本) + text_specs = {t[0]: t for t in TEXT_TRAITS} + ssc_p = params["ssc"] + obs_rows = [] + for info in combo_info: + r = info["rounds"].get("line") + if not r or not r.get("active") or "sel" not in r: + continue + obs_year = r["year"] + op = personnel_ids[obs_year % len(personnel_ids)] + bv = info["bv"] + for ei in r["sel"]: + ei = int(ei) + tid = info["tree_ids"][ei] + ssc_val = clip_round(ssc_p["mu"] + bv[ei, ssc_p["idx"]] + pyrand.gauss(0, ssc_p["sd"] * 0.4), + ssc_p["vmin"], ssc_p["vmax"]) + obs_rows.append({**_base(), + "tree_id": tid, "trait_id": trait_id.get("ssc"), "obs_date": date(obs_year, 7, 20), + "obs_year": obs_year, "obs_value": f"{ssc_val:.1f}", "obs_type": "numeric", + "operator_id": op, "status": 1, "remark": "模拟数据"}) + for tcode, t_tid in (("flesh_color", trait_id.get("flesh_color")), + ("fruit_shape", trait_id.get("fruit_shape"))): + val = pyrand.choices(text_specs[tcode][4], weights=text_specs[tcode][5], k=1)[0] + obs_rows.append({**_base(), + "tree_id": tid, "trait_id": t_tid, "obs_date": date(obs_year, 7, 20), + "obs_year": obs_year, "obs_value": val, "obs_type": "text", + "operator_id": op, "status": 1, "remark": "模拟数据"}) + if obs_rows: + counts["observations"] += await insert_chunked(session, ObservationModel, obs_rows) + await session.commit() + + # 5) 试验处理:每个研究点 2-3 条 + study_ids = [r[0] for r in (await session.execute(select(TrialStudyModel.id))).all()] + treat_rows = [] + for sid in study_ids: + for tf, lv, desc in (("栽培管理", "对照", "常规栽培"), ("栽培管理", "提质", "控产提质"), + ("施肥", "常规", "常规施肥"), ("施肥", "减量", "减量增效")): + if pyrand.random() < 0.65: + treat_rows.append({**_base(), "trial_study_id": sid, "factor": tf, + "level": lv, "description": desc, "remark": "模拟数据"}) + if treat_rows: + await insert_chunked(session, TreatmentModel, treat_rows) + counts["treatments"] += len(treat_rows) + + +# --------------------------------------------------------------------------- +# 汇总与自检 +# --------------------------------------------------------------------------- +async def report_counts(session, counts) -> None: + print("[sim] 生成统计:") + print(f" 杂交组合 : {counts['combos']}") + print(f" 单株 : {counts['trees']}") + print(f" 单株评价 : {counts['evals']}") + print(f" 性状观测 : {counts['obs']}") + print(f" 农事操作 : {counts['field_ops']}") + print(f" 选育结果(晋级+淘汰) : {counts['selections']}") + print(f" 定植批次 : {counts['plantings']}") + print(f" 区域试验/研究点/参试: {counts['trials']}/{counts['studies']}/{counts['entries']}") + for t, key in (("bre_clone", "clones"), ("bre_germplasm", "germplasm"), ("bre_pedigree", "pedigree")): + n = (await session.execute(text(f"SELECT count(*) FROM {t}"))).scalar() + counts[key] = n + print(f" {t:<16} : {n}") + for t, key in (("bre_pollination", "pollination"), ("bre_pollen", "pollen"), + ("bre_seed_lot", "seed_lots"), ("bre_seed_treatment", "seed_treatments"), + ("bre_seedling", "seedlings"), ("bre_propagation", "propagations"), + ("bre_environment_condition", "env_conds"), ("bre_observation", "observations"), + ("bre_treatment", "treatments")): + print(f" {t:<20} : {counts[key]}") + + +async def verify(session, counts, dry: bool = False) -> None: + print("[sim] 自检 ...") + problems = [] + min_obs = 8000 if dry else 300000 + + def chk(name, cond): + ok_ = bool(cond) + print(f" [{'OK ' if ok_ else 'BAD'}] {name}") + if not ok_: + problems.append(name) + + n_tree = (await session.execute(text("SELECT count(*) FROM bre_tree"))).scalar() + n_obs = (await session.execute(text("SELECT count(*) FROM bre_trait_observation"))).scalar() + n_obs_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_trait_observation o LEFT JOIN bre_tree t ON t.id=o.tree_id " + "WHERE o.tree_id IS NOT NULL AND t.id IS NULL"))).scalar() + n_eval_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_tree_evaluation e LEFT JOIN bre_tree t ON t.id=e.tree_id " + "WHERE t.id IS NULL"))).scalar() + n_sel_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_selection_result s LEFT JOIN bre_tree t ON t.id=s.tree_id " + "WHERE t.id IS NULL"))).scalar() + n_dup_tree = (await session.execute(text( + "SELECT count(*) FROM (SELECT tree_no FROM bre_tree GROUP BY tree_no HAVING count(*)>1) x"))).scalar() + n_dup_clone = (await session.execute(text( + "SELECT count(*) FROM (SELECT clone_code FROM bre_clone GROUP BY clone_code HAVING count(*)>1) x"))).scalar() + n_germ_dup = (await session.execute(text( + "SELECT count(*) FROM (SELECT cultivar_name FROM bre_germplasm " + "WHERE cultivar_name ~ '\\-\\d{4}$' GROUP BY cultivar_name HAVING count(*)>1) x"))).scalar() + n_sp = (await session.execute(text("SELECT count(*) FROM bre_tree WHERE stage='sp'"))).scalar() + n_ap = (await session.execute(text("SELECT count(*) FROM bre_tree WHERE stage='ap'"))).scalar() + n_line = (await session.execute(text("SELECT count(*) FROM bre_tree WHERE stage='line'"))).scalar() + n_elim = (await session.execute(text("SELECT count(*) FROM bre_tree WHERE status='eliminated'"))).scalar() + n_sel_tree = (await session.execute(text( + "SELECT count(*) FROM bre_tree WHERE status IN ('selected','primary','key','preserved')"))).scalar() + n_line_germ = (await session.execute(text( + "SELECT count(*) FROM bre_germplasm WHERE cultivar_name ~ '\\-\\d{4}$'"))).scalar() + n_sel_total = (await session.execute(text("SELECT count(*) FROM bre_selection_result"))).scalar() + n_polli = (await session.execute(text("SELECT count(*) FROM bre_pollination"))).scalar() + n_pollen = (await session.execute(text("SELECT count(*) FROM bre_pollen"))).scalar() + n_lot = (await session.execute(text("SELECT count(*) FROM bre_seed_lot"))).scalar() + n_treat = (await session.execute(text("SELECT count(*) FROM bre_seed_treatment"))).scalar() + n_seedl = (await session.execute(text("SELECT count(*) FROM bre_seedling"))).scalar() + n_prop = (await session.execute(text("SELECT count(*) FROM bre_propagation"))).scalar() + n_env = (await session.execute(text("SELECT count(*) FROM bre_environment_condition"))).scalar() + n_obs_g = (await session.execute(text("SELECT count(*) FROM bre_observation"))).scalar() + n_trt = (await session.execute(text("SELECT count(*) FROM bre_treatment"))).scalar() + n_seedl_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_seedling s LEFT JOIN bre_seed_lot l ON l.id = s.seed_lot_id " + "WHERE s.seed_lot_id IS NOT NULL AND l.id IS NULL"))).scalar() + n_prop_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_propagation p LEFT JOIN bre_clone c ON c.id = p.produced_clone_id " + "WHERE p.produced_clone_id IS NOT NULL AND c.id IS NULL"))).scalar() + n_obs_tree_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_observation o LEFT JOIN bre_tree t ON t.id = o.tree_id " + "WHERE o.tree_id IS NOT NULL AND t.id IS NULL"))).scalar() + n_clone_single = (await session.execute(text( + "SELECT count(*) FROM (SELECT clone_id FROM bre_tree " + "WHERE clone_id IS NOT NULL AND is_deleted = false " + "GROUP BY clone_id HAVING count(*) < 2) x"))).scalar() + n_ramet_orphan = (await session.execute(text( + "SELECT count(*) FROM bre_tree t LEFT JOIN bre_clone c ON c.id = t.clone_id " + "WHERE t.clone_id IS NOT NULL AND c.id IS NULL"))).scalar() + + chk("单株数 > 0", n_tree > 0) + chk(f"观测数 ≥ {min_obs / 10000:g}万", n_obs >= min_obs) + chk("观测无孤儿 tree_id", n_obs_orphan == 0) + chk("评价无孤儿 tree_id", n_eval_orphan == 0) + chk("选育结果无孤儿 tree_id", n_sel_orphan == 0) + chk("单株编号无重复", n_dup_tree == 0) + chk("克隆编号无重复", n_dup_clone == 0) + chk("晋级种质无重名", n_germ_dup == 0) + chk("晋级阶段分布合理", n_sp >= n_ap >= n_line and n_sp > 0) + chk("有淘汰树且选中树存活", n_elim > 0 and n_sel_tree > 0) + chk("晋级种质(line+)> 0", n_line_germ > 0) + chk("选育结果 = 晋级 + 淘汰", n_sel_total == counts["selections"]) + chk("授粉/花粉/种子链有数据", n_polli > 0 and n_pollen > 0 and n_lot > 0 and n_treat > 0 and n_seedl > 0) + chk("扩繁/气象/通用观测/处理有数据", n_prop > 0 and n_env > 0 and n_obs_g > 0 and n_trt > 0) + chk("育苗无孤儿 seed_lot_id", n_seedl_orphan == 0) + chk("扩繁无孤儿 clone_id", n_prop_orphan == 0) + chk("通用观测无孤儿 tree_id", n_obs_tree_orphan == 0) + chk("每个克隆 ≥ 2 棵树(ABLUP G1 可建模)", n_clone_single == 0) + chk("克隆树无孤儿 clone_id", n_ramet_orphan == 0) + + print(f" [信息] stage 分布: seedling/sp/ap/line = " + f"{n_tree - n_sp - n_ap - n_line}/{n_sp}/{n_ap}/{n_line}") + if problems: + print(f"[sim] 自检失败 {len(problems)} 项: {problems}") + raise SystemExit(1) + print("[sim] 自检全部通过 ✓") + + +if __name__ == "__main__": + asyncio.run(main()) diff --git a/doc/桃育种系统业务链路与数据模型说明书 v1.1.docx b/doc/桃育种系统业务链路与数据模型说明书 v1.1.docx new file mode 100644 index 0000000..4643203 Binary files /dev/null and b/doc/桃育种系统业务链路与数据模型说明书 v1.1.docx differ diff --git a/doc/桃育种系统统计引擎实施记录.md b/doc/桃育种系统统计引擎实施记录.md index a599074..3416610 100644 --- a/doc/桃育种系统统计引擎实施记录.md +++ b/doc/桃育种系统统计引擎实施记录.md @@ -1,8 +1,8 @@ # 桃育种系统 · 统计引擎实施记录 > **版本**:v1.4(2026-08-04 追加统计严谨·预测完整性与 MLOps 复现性:PA 落库 / 溯源 / germplasm_id 填充 / 砧木扩列 / Calo 相关校正 / 空间竞争协变量 + 组合得失漏斗) -> **定位**:记录统计引擎多轮 P0/P1 落地 + 田间刚需(花粉档案)+ 决策正确性(无性系级指数)+ 统计严谨(配合力交配设计 / 相关校正 / 空间竞争)+ MLOps 复现性 + 组合得失漏斗的实现、验证与文件清单——**G×E 交互引擎**(2026-08-03)、**性状方向标注**(2026-08-04)、**花粉档案**(2026-08-04)、**无性系级指数**(2026-08-04)、**配合力交配设计**(2026-08-04)、**统计严谨·MLOps 复现性**(2026-08-04,①Calo②PA③溯源④germplasm_id⑤砧木⑥空间竞争)与**组合得失漏斗**(2026-08-04,⑧ v_combination_funnel)。本文档为持续追加的落地台账,与规格文档 `桃育种系统模块扩展需求规格.md`(v2.10)配套;规格管"应该长什么样",本文管"实际落地了什么、怎么验证的"。 -> **相关探针**:`Temp/claude/e2e_gxe_20260803.py`(G×E,5 场景)· `Temp/claude/e2e_direction_20260804.py`(方向标注,7 组)· `Temp/claude/e2e_pollen_20260804.py`(花粉档案)· `Temp/claude/e2e_clone_index_20260804.py`(无性系级指数,6 组)· `Temp/claude/e2e_combining_design_20260804.py`(配合力交配设计,5 组)· `Temp/claude/e2e_prediction_rigor_20260804.py`(②③④,4 组)· `Temp/claude/e2e_rootstock_20260804.py`(⑤)· `Temp/claude/e2e_corr_spatial_20260804.py`(①⑥,3 组)· `Temp/claude/e2e_combination_funnel_20260804.py`(组合得失漏斗,3 组)。 +> **定位**:记录统计引擎多轮 P0/P1 落地 + 田间刚需(花粉档案)+ 决策正确性(无性系级指数)+ 统计严谨(配合力交配设计 / 相关校正 / 空间竞争)+ MLOps 复现性 + 组合得失漏斗的实现、验证与文件清单——**G×E 交互引擎**(2026-08-03)、**性状方向标注**(2026-08-04)、**花粉档案**(2026-08-04)、**无性系级指数**(2026-08-04)、**配合力交配设计**(2026-08-04)、**统计严谨·MLOps 复现性**(2026-08-04,①Calo②PA③溯源④germplasm_id⑤砧木⑥空间竞争)与**组合得失漏斗**(2026-08-04,⑧ v_combination_funnel)。本文档为持续追加的落地台账,与规格文档 `桃育种系统模块扩展需求规格.md`(v2.10)配套;规格管"应该长什么样",本文管"实际落地了什么、怎么验证的"。 配套《桃育种系统业务链路与数据模型说明书》v1.1(2026-08-06 重生成,已同步 O1 花期数据消费补强 / N4 plan→trial→target 血缘 / 亲本校验 / 世代闭环落库 / 新增 §九 QA 与 CI 章节)与代码保持一致,三者(规格 + 实施记录 + 说明书)互为索引。 +> **相关探针**:`Temp/claude/e2e_gxe_20260803.py`(G×E,5 场景)· `Temp/claude/e2e_direction_20260804.py`(方向标注,7 组)· `Temp/claude/e2e_pollen_20260804.py`(花粉档案)· `Temp/claude/e2e_clone_index_20260804.py`(无性系级指数,6 组)· `Temp/claude/e2e_combining_design_20260804.py`(配合力交配设计,5 组)· `Temp/claude/e2e_prediction_rigor_20260804.py`(②③④,4 组)· `Temp/claude/e2e_rootstock_20260804.py`(⑤)· `Temp/claude/e2e_corr_spatial_20260804.py`(①⑥,3 组)· `Temp/claude/e2e_combination_funnel_20260804.py`(组合得失漏斗,3 组)。 · `Temp/claude/e2e_o1_bloom_20260805.py`(O1 花期数据消费补强:四态软警示 + 窗口覆盖三态 + 库存五态 + 反向门禁,4 组) --- diff --git a/doc/桃育种系统规划对照总表.md b/doc/桃育种系统规划对照总表.md index 428a94b..114dece 100644 --- a/doc/桃育种系统规划对照总表.md +++ b/doc/桃育种系统规划对照总表.md @@ -2,7 +2,7 @@ > **版本**:2026-07-29 新增 > **定位**:统一三份文档(业务功能规划报告 / 分层路线图 v1.1 / 模块扩展需求规格 v2.11)的模块级映射,消除口头"指回§3/§4"。 -> **状态约定**:✅ 已落地(基础字段) | ⏳ 待建(规划已定,待生成器落地) | ⏳ 待灌(数据) | — 引擎/端点(不建表)。**2026-08-04 追加**:① 统计引擎两轮 P0 已代码落地(G×E + 性状方向标注);② 花粉档案 §3.19 落地(bre_pollen + 授粉窗口校验);③ 选择指数升级为**无性系级聚合**(决策正确性,可靠加权聚合 + 自株退化 + 系级决选);④ 配合力 GCA/SCA 支持**交配设计 `design_type`**(统计严谨:完全双列/部分双列/line×tester(NCII)/NCIII,`combining.solve` 按设计分支双因素模型 + 运行按设计过滤组合);⑤ **统计严谨·MLOps 复现性**(① Smith-Hazel 遗传相关改 **Calo 可靠性校正** r_g=r_EBV/√(rel_i·rel_j);② **PA 落库** pa=√reliability + 批次 accuracy 真 PA;③ **MLOps 溯源** data_version/input_hash/engine_version(同数据重跑哈希一致、改观测必变);④ **germplasm_id 填充**(亲本级 EBV 查询可用);⑤ **砧木字典扩列** 矮化类/砧穗亲和性;⑥ **空间竞争协变量** covariate=competition(Chebyshev 邻域株数,无坐标显式报错));⑥ **组合得失漏斗 v_combination_funnel**(按组合汇总 花→果→种→苗→定植→树→入选 六级派生指标——结实率/出苗率/选择强度滚到组合层,配合力第一手证据,只读端点可查);⑦ **分子育种「建数据能力」**(§3.7 四表落地:bre_marker / bre_genotyping_dataset / bre_genotype_sample(weld 补 sample_name)/ bre_genotype_call,全量 CRUD + Excel 导入导出 + **VCF 导入**(marker/sample 自动建号、GT 编码、幂等)+ **DNA 指纹防错**(无性系混杂检测,同源不一致/异源一致两类标记));⑧ **精修项四连发**(① **DUS/品种保护「数据能力(三表)」** bre_dus_descriptor/test/observation + UPOV TG/53 桃描述符种子 15 条 + 菜单 900230 真实页与 G 段按钮 900800-900808;② **stage 感知** bre_trait.stage(童期/成株)驱动 BLUP/选择指数/决策预览的过滤与警示 + 批次 stage 落库;③ **k-fold 交叉验证外部验证** ABLUP+GXEBLUP 双支持(sire 家系分层留出、逐折容错,bre_cv_result/fold 落库,MLOps 溯源复用);④ **按 site EBV 公平性报告** site 系统性偏差 deviation/flagged + 排名一致性 Spearman + G×E 模式);⑨ **六项能力完善**(① **MT-BLUP 遗传相关**——成对双性状 REML(`mtblup.solve_bivariate`:y 堆叠 Var(u)=G0⊗A,Va/Ve 取单性状 REML、仅对 ρ 黄金分割求极大精确 REML),`run_genetic_corr` 逐对 bivariate 组装 G0 + Higham 半正定投影,落 `bre_genetic_corr_result`(job_type=GENCORR),Smith-Hazel 增 `g_method="mtblup"`(单对不收敛回退 Calo+warning、G 非正定特征值截断、`g_source` 落库可审计);② **GBLUP/ssGBLUP 引擎**——`genomic.py`:VanRaden G method1/2 + MAF 过滤 + blend 岭;`solve_gblup` 仅基因型株、`solve_ssgblup` 单步法 H⁻¹=A⁻¹+[[0,0],[0,G⁻¹−A22⁻¹]] 且 **V 块用 H 协方差**;`run_gblup` 样本 tree 直连/名称兜底、未映射+多等位跳过+warning,落 PredictionModel(GBLUP/ssGBLUP)+EBV 行;③ **DUS 特异性统计检验**——QN LSD(t_crit 由 `fdist.f_icdf` 二分反解 F(1,df))+ PQ/QL 状态众数、**必测描述符**驱动建议、参照=同 trial_study 自动 + 显式覆盖、落 `analysis_json` + conclusion 自动建议;④ **AMMI/Finlay-Wilkinson 稳定性**——FW 环境指数回归 b/R²/se_b/flag_stable(完美拟合 b≈1 判稳)、AMMI 残差 SVD→IPC/ASV/ecovalence、`run_stability` site/year 两维、格子 <2×2→409,落 `bre_stability_result`;⑤ **MLOps 重训+版本回滚**——`bre_prediction.is_active` 版本链(无 active 时首批 active)、`model/drift` 输入哈希漂移检测(轻量 `_gather_digest_inputs`)、`model/retrain` 漂移才重训+active 转移、`model/activate` 回滚原语;⑥ **BLUP stage 拆分**——`bre_trait_observation.stage` 观测级发育阶段、`run_ablup(stage)` 按 coalesce(obs.stage,trait.stage) 过滤+model_name `_{stage}` 后缀+批次 stage 落库、未指定但有分歧→note 警示);⑩ **田间试验精度补强·砧木×接穗随机互作(G×R)**——rootstock 由固定哑变量升级为**第三随机效应槽**(`run_ablup(gxr=True)` → method=GXRBLUP/job_type=GXR,按 `(基因型,砧木)` 分组建 Z₂ 随机互作、克隆坍缩系谱与 G×E 共享、rootstock 强制移出固定效应防共线、不可辨识门禁(同基因型无跨砧木/单元内无重复/残差 df<1)→409、σ²gxr/gxr_ratio/n_cross_rootstock 落批次 note;**Z₂ 单槽位 → G×E/G×R 互斥**(双开 409 + 前端互斥 watcher);求解器 `_solve_gxe` 加 factor_label/factor_name(G×E 警告文案逐字节不变、G×R 显「砧木」)、ENGINE_VERSION 1.2.0;G×E tab 加 G×R 复选框、批次表同显 σ²gxr);⑪ **基因组选择(GS)实证严谨性·GBLUP/ssGBLUP 专用 k-fold 交叉验证**——`run_cv` 加 `dataset_id/method/maf_min`(method=KFCV/GBLUP·KFCV/ssGBLUP),`genomic.kfold_cv_genomic` 掩蔽留出 GEBV vs 表型实证 pearson/RMSE(G/H/Hinv 只建一次、固定种子随机分层——GS 同世代样本无家系树,区别于 ABLUP 的 sire 家系留出)、`_gather_gblup_inputs` 提取(run_gblup↔GS CV 复用)、`_build_h` 缺失基因型补 base、ENGINE_VERSION 1.1.0、_DATA_VERSION v2.11;SSR 多等位 dummy 编码留待后续);⑫ **现代桃育种领域覆盖·S-等位基因交配兼容性 + 抗病/需冷量性状字典**(① `bre_germplasm.s_alleles`(String32,Sf=自交亲和型)+ `cross_combination` create/update 兼容校验——任亲本含 Sf→放行、共享 2 个 S 等位→409「配了不结」硬门禁、共享 1 个→`out.s_compat` 半兼容警示;②③ `bre_trait` 种子 5 条 numeric(细菌性穿孔病/褐腐病/白粉病 0-5 级病情指数 + 需冷量 h + 需热量 GDD,category「抗病性」「生态适应性」,into_ebv=1、direction=asc 进选种目标——统计引擎仅消费 numeric);④ **MAS/QTL/GWAS 闭环**(GWAS 引擎 GLM+PC:SVD 群体结构校正 + 逐标记单标记回归 + Bonferroni/BH-FDR + QTL 区间合并,纯 numpy 复用 fdist 零新依赖;bre_gwas_result/bre_gwas_snp/bre_qtl/bre_mas_panel/bre_mas_panel_marker 五新表;服务端 gwas/run|list|{id} + qtl CRUD + mas-panel CRUD(选标记 favorable_dose/direction/effect);decision_preview 加 marker 条件——童期幼苗无表型/无 EBV 也可被标记辅助选入;前端 900214 GWAS/QTL 页翻转(Manhattan/QQ 裸 echarts + 运行对话框 + 三 tab);e2e 7 段全过 + GBLUP/ssGBLUP 回归零回归,2026-08-04);⑭ **九缺口补齐(2026-08-04,用户重扫 12 项缺口剔除 A1/A4 余 9 项拍板全做,三批独立 e2e + 回归)**——批1 经典侧计算端点(纯计算不建表):**数据质量** `POST /statistics/data-quality`(缺失率/CV + IQR + MAD 稳健 z-score 双法异常标记)+ **遗传增益** `POST /statistics/genetic-gain`(ΔG=k·r_g·σ_A,k 经 Acklam 有理近似 Φ⁻¹ 纯 numpy 零 scipy)+ **主动选配** `POST /statistics/mating-recommend`(S-等位硬过滤 + 亲缘 A 惩罚 score=w_ebv·mid_parent_EBV−w_kin·max(0,r−threshold),不落库);批2 分子侧严谨性:**EMMAX 混合模型 GWAS**(`gwas._emmax`:零模型方差分量复用 `genomic._profile_solve` + eigh 对角化 K→V⁻¹ + 逐标记 GLS t²→fdist,method=emmax,run_gwas 后处理抽共享函数 GLM 逐字节不变)+ **SSR 多等位虚拟编码**(VCF 按 ALT 数定 ssr/snp、`_allelic_from_gt` 每等位一列、`build_g_matrix` 多等位 VanRaden 推广)+ **QTL×E 环境互作**(`POST /statistics/gwas-qtl-x-e` 分环境 GWAS 合并 stable/env-specific/异号 G×QTL 警示,不建表)+ **MAS 面板语义**(bre_mas_panel_marker.mode/favorable_allele/haplotype_group 三列 + `_mas_hit_map` additive/dominance/recessive/allele/haplotype 五语义);批3 求解器核心:**AR1×AR1 空间协方差**(`blup.solve_spatial`:e~N(0,σ²e·R)、R=AR1(ρ)⊗AR1(ρ) R⁻¹ 纯 numpy eigh、REML 对 (h²,ρ) 坐标上升黄金搜索、run_ablup(spatial=True)→method=AR1×AR1、与 G×E/G×R 互斥、缺坐标 409)+ **稀疏 A⁻¹+共轭梯度**(`_build_ainv_sparse` COO + `_cg_solve`、n>N_SPARSE(1000) 阈值分派 solver=sparse-cg、Hutchinson k=100 估 PEV 对角 reliability 近似、稠密路径逐字节不变)。`_DATA_VERSION` v2.13、blup 1.3.0 / genomic 1.2.0 / gwas 1.1.0。详见 `桃育种系统统计引擎实施记录.md`);⑮ **十二项真缺口补齐(2026-08-05,用户全量现状盘点收敛 12 件真缺口、剔除过时「待建」标记,拍板只补这 12 件,分三批独立 e2e + 全量回归零回归)**——批1 引擎三件(solver 扩展独立路径):**ssGWAS**(`gwas._ssgwas`:`build_g_matrix` 岭保正定构 G → `solve_gblup` 全样本 GEBV û → **Wang et al. 单步 GWAS MEM 对角近似**反推标记效应 `ê_j=(M_j'Z'G⁻¹û)/(2Σp_j(1−p_j))`、Wald t=ê/SE → fdist,共享 `_finalize`/`_bh_qvalues`/`_cluster` 后处理,method=ssgwas 透传+前端下拉,gwas 1.2.0)+ **全 MT-BLUP solve_multi**(m 性状 y 堆叠 Var(u)=G0⊗A、**EM-REML 迭代** G0+ve 单性状初始+特征值截断/Higham 投影保正定、max_iter 兜底,输出 G0 全元素+r_g 矩阵+n_iter/converged/warning,`run_genetic_corr(full_mtblup=True)` 一次 m 性状给完整 G0,逐对 bivariate 保活,mtblup 1.1.0)+ **AR1×AR1 双参 ρ aniso**(R_ij=ρ_row^|Δrow|·ρ_col^|Δcol|、ρ_row/ρ_col 交替坐标上升,`run_ablup(spatial_aniso=True)` → method=AR1×AR1(aniso)、note 落 rho_row/rho_col,**aniso=False 与 v1 逐字节一致**,blup 1.4.0);批2 模块四件(标准 CRUD 五件套+前端+菜单翻转+文案清理):**field_operation 农事操作 900051**(op_type 字典施肥/喷药/修剪/灌溉/疏果/套袋/采收 bre_dict 灌入)+ **observation 通用观测 900050**(EAV plot 级通用录入、obs_type numeric/text/date)+ **breeding_report 育种报告 900075**(CRUD + `POST /report/generate` 聚合 bre_prediction/bre_stability_result/bre_cv_result 结构化报告)+ **analysis_dataset 分析数据集 900070**(新表 bre_analysis_dataset jsonb×3 + status,「查看性状值」跳 statistics trait-values 复用既有视图);批3 底座五件:**审计切面**(base_crud 三写注入 CREATE/UPDATE 字段级 diff/DELETE + 批量导入仅汇总 IMPORT 一条 + `GET /bre/audit/list`)+ **统一编号服务**(`app/utils/number_gen.py` NumberGenService + `bre_sequence` + **pg_advisory_xact_lock 原子取号**,迁移 cross_combination/tree/selection_result/seedling 4 处格式保持)+ **超期预警定时任务**(每日 job 按 selection_rule 扫描树 N 年未决策/花粉批次过期/seed_lot 超期 → `bre_alert` + `GET /bre/alert/list` + 前端「预警中心」)+ **备份容灾**(`backend/scripts/pg_backup.sh`/`pg_restore.sh`,Temp/backup 按日期 + N 天轮转,每日 02:00 系统 job)+ **BrAPI + MIAPPE 互操作**(`/brapi/v2` 只读五端点 germplasm/observationvariables/studies/observations/miappe,BrAPI 2.x metadata+pagination+result 包裹,bre_germplasm→germplasmDbId/genus=Prunus/species=persica/subtaxa 映射);`_DATA_VERSION` v2.14、gwas 1.2.0 / mtblup 1.1.0 / blup 1.4.0。详见规格 §8.24–8.26 / `桃育种系统统计引擎实施记录.md` §十七);⑯ **桃育种业务流全面完善(2026-08-05,用户要求「对现有功能,结合桃育种实际业务流程与需求,进行全面测试分析,并完善」,三路审计收敛 10 项真缺口,拍板全部实施 + 两裁决:**导入按钮维持隐藏**(23+ 页显式隐藏非缺口)/ **采收入口走农事操作**(op_type=harvest 结构化产量表单,非独立模块),分三批独立 e2e + 全量回归零回归)**——批A 业务数据链(A1-A6):**果实采收结构化入口**(bre_field_operation 加 yield_kg/fruit_count/avg_fruit_weight/marketable_rate,op_type=harvest 树级采收同步写 bre_trait_observation 4 行、plot 级只留操作记录)+ **产量构成/物候/品质/生长量/描述性状种子 19 条**(bre_yield_phenology_traits.sql:单株产量/果数/均果重/好果率/裂果率 + 盛花期/果实成熟期/落叶期 + 可溶性固形物/可滴定酸/固酸比 + 主干周长/树高/冠幅 + 果形/茸毛/离核/果皮底色/着色类型)+ **种质级观测**(bre_observation.germplasm_id 描述数据入口)+ **授粉→收种→种子处理补链**(seed_lot.pollination_id + seed_treatment.seed_lot_id)+ **EAV 观测校验**(numeric/date parse + valid_min/max 越界 + 同(tree/plot/germplasm,trait,year,date)重复 409)+ **育种阶段自动设置 + 成株性状门禁**(tree stage 按 breeding_generation 推断、童期树录成株性状 409);批B 统计算法(B1-B4):**Type-B 多环境遗传力**(run_type_b_heredity 复用 mtblup.solve_bivariate 把环境当"性状"、env_dim=site/year、落 bre_type_b_result + /type-b-heredity + /type-b/{id})+ **遗传相关 genetic 模式**(correlation mode=genetic 逐对 bivariate 组装 G0,pheno 分支逐字节不变、不落库)+ **选择指数自动权重**(auto_weights=True → default_h2 兜底 0.1、强制 use_h2=False、direction 翻转)+ **UPGMA 层次聚类**(run_cluster 纯 numpy,POST /statistics/cluster 不落库);批C 前端(C1-C4):观测录入增强(observation plot/站点/germplasm 过滤 + germplasm 表单、field_operation harvest 结构化表单)+ 统计可视化 ECharts(correlation 热图 + gblup 可靠性热图/趋势折线 + type-B tab + cluster 树状图)+ combination-funnel 前端入口对话框 + stage_phenotype 容器同步(通用观测/农事操作翻真实页);`_DATA_VERSION` v2.15。详见规格 §8.27 / `桃育种系统统计引擎实施记录.md` §十八);⑰ **遗传链完整性两次修复(2026-08-05,用户先后报告「建种质断链」与「近交惩罚静默关闭」,核实确认 + 拍板落地,e2e 验证通过)**——① **世代闭环修复**:`selection_result._promote_tree_to_clone` 单株晋级到 line/regional_trial/released 建种质时**同步落 `bre_pedigree`**(child_code=品种名、dam/sire=tree.dam/sire 回退组合 female/male_parent_id、combination_id、generation、child_code 幂等)→ 晋升种质不再被统计引擎当 founder,跨世代系谱链闭合(判别点 A[T2,T3]=0→0.125、A[gT1,gA]=0.5);② **近交惩罚「静默关闭」修复**:`mating_recommend` 亲缘矩阵旧靠 `bre_pedigree.child_code` 精确字符串匹配建系谱,名称大小写/空格/代号错配即全空 → 所有配对 r=0 → kin_pen=0 → 近交规避失效且无告警;重写 `_germplasm_pedigree` **多源级联**(① 规范化 strip+casefold 名称反查显式系谱为权威 → ② `tree.germplasm_id` 直连 FK 兜底(多对去重后唯一一致才采纳、歧义回退 founder)→ ③ founder;祖先 BFS 展开 + 去重,修复候选既是 child 又是祖先时系谱重复抛『系谱个体重复』崩溃);`mating_recommend` 新增 **kinship_status(ok/partial/none) + kinship_warning + 逐对『亲缘未解析(r 按 0 计)』旗标**,绝不静默 no-op;设计红线 **FK 只兜底不覆盖显式系谱**(A6 回归 r=0.5 保真)。e2e 14 组断言 + 回归 A6/batch/gaps_base 零回归 + HTTP 核验新字段。详见规格 §8.28 / `桃育种系统统计引擎实施记录.md` §18.5/§18.6。);⑱ **选配「亲本 EBV」回退逻辑方法学修复(2026-08-05,用户报告 `mating_recommend` 旧实现对缺直接 EBV 的候选用「其子代树 EBV 均值」代表亲本育种值——① progeny-test 值非亲本自身值;② 双亲都缺时同一子代树等权计入 dam 与 sire 两侧稀释;③ 跨组合混合无区分,会误导选配排序。核实确认 + 拍板落地)**——`mating_recommend` EBV 回退重写:**direct(germplasm_id 直连亲本自身预测值)优先**、有 direct 绝不被子代树值污染;回退仅当无自身值时用、**显式标注「子代测验近似」并降权 0.5 参与打分**(`_EBV_SOURCE_FACTOR={direct:1.0,progeny:0.5,missing:0.0}`,mid=(f_a·ebv_a+f_b·ebv_b)/2,回退不掩盖 direct);**双亲均缺直接值的子代树剔除**(防双侧等权重复计);记 n_progeny/n_combos 辨别跨组合;返回体加 `candidate_ebv` + `ebv_source.coverage` + 逐对旗标「亲本EBV为子代测验近似(降权0.5)」/「亲本无EBV(按0计)」——EBV 来源显式可查绝不静默降级。e2e 5 组断言 + 回归 A6/亲缘/batch/gaps_base 零回归 + HTTP 核验(真实数据 coverage={direct:0,progeny:0,missing:12})。详见规格 §8.29 / `桃育种系统统计引擎实施记录.md` §18.7。);⑲ **数据质量报告「离群值」方向语义错位修复(2026-08-05,用户报告 §8.23 数据质量端点 IQR/MAD 把任何偏离分布的值都标「离群株」,与育种目标方向错位——正向偏离正是要选的精英单株,核实确认 + 拍板落地)**——`data_quality_report` 方向感知分流:按 `bre_trait.direction`(desc/None=高值优、asc=低值优)判定 `desirable_high`,IQR/MAD 命中的株按偏离侧**双通道分类**(`is_high=v>med`、`is_elite=is_high==desirable_high`):**同向偏离 → `extreme_candidates` 极值候选(精英)**(class=extreme_candidate)/ **反向偏离 → `outliers` 疑似录入错误**(class=outlier,高值优的异常低值、低值优的异常高值),两表互斥;返回体加 `direction` + `extreme_candidates` + `summary.n_extreme_candidates`、`has_warning/n_outliers` 语义收紧为疑似录入错误计数;前端数据质量 tab 改双表(精英 success / 错误 danger)+ 选育方向指示 + 无误报绿色 info alert。e2e 17 株断言(精英 16.0 进 extreme_candidates/side=high、错误 2.0 进 outliers/side=low、互斥不串表)+ 回归 A6/亲缘/batch/gaps_base 零回归 + HTTP 核验(asc 裂果率 bB_DIR5668 与 desc 果重 dA_DIR5668 方向/双表/计数经真实端点正确返回,真实数据无极端值→双表空无误报)。详见规格 §8.30 / `桃育种系统统计引擎实施记录.md` §18.8。);⑳ **ΔG 投影假设声明 + 工程健壮性打磨(2026-08-05,用户报告两点均非功能 bug,核实确认 + 拍板落地,e2e 验证通过)**——① **ΔG 投影假设显式声明**:`genetic_gain` 旧实现 `next_pheno_mean = pheno_mean ± dg` 隐含「环境均值恒定」,年际气候漂移会使投影偏;属投影演示性质,返回体新增 `note` 假设声明(环境均值恒定、未计入年际气候漂移与栽培措施变化;ΔG=k·r_g·σ_A 仅含遗传增益分量,next_pheno_mean 为投影演示值而非实测预测),前端 gain tab info alert 展示;② **空间模型缺坐标优雅降级**:`run_ablup(spatial=True)` 旧实现任一观测树缺 row/col 即整批 ValueError,现改为**剔除缺坐标株表型(保留系谱 → EBV 仍由系谱预测)、仅全缺坐标才整批 409**,剔除数显式写批次 note「缺坐标剔除N株表型(EBV由系谱预测)」、input_hash 用剔除后子集;③ **numpy 工程清理**:7 处函数体内重复 `import numpy` 收敛为模块级单次 + `data_quality_report` 逐株 flag 循环改 numpy 布尔掩码向量化(输出逐字节不变)。e2e spatial 降级场景(5 株含 1 缺坐标 → method=AR1×AR1/train_n=4/note 含剔除/5 株均有 EBV)+ A6 note 断言 + 回归 spatial/A6/亲缘/batch/gaps_base 五套零回归 + vue-tsc EXIT=0 + HTTP 核验(genetic-gain 返回 note、data-quality 返回 direction)。详见规格 §8.31 / `桃育种系统统计引擎实施记录.md` §18.9。);㉑ **A–F 六项功能级缺口·批1 父本验证 + 近交衰退(2026-08-05,用户全模块 grep 实证 A–F 六项缺口均为「有无/深度」问题、非正确性 bug,逐项复核确认 + 拍板「全部六项分批落地」;批1 独立 e2e + 全量回归零回归)**——**A 父本验证(consistency 不落库)**:① `pollination` 花粉来源一致性——`_validate_pollen_parent`(create/update 统一):花粉批次声明采集父本树(pollen.male_tree_id)时授粉声明父本树须一致否则 409(含「不一致」文案),批次无采集父本(混合/外地)或授粉未声明放行不武断;② `tree` 显式亲本 vs 组合父母本警示——`_resolve_parents` 返回 `(dam,sire,warning)`:显式 dam/sire 偏离组合 female/male_parent 时返回体 `warning`「显式母本/父本与组合母本/父本不一致(以显式为准)」(双侧/单侧聚合),**warning 不阻断**(tree 可独立建,武断 409 误伤),schema 只读字段不落库;③ `cross_combination` 同质自交门禁——`_check_parent_roles` 中 `female_parent_id==male_parent_id` → 409「母本与父本同为「{名}」,自交组合请确认」(create/update 生效亲本统一,与 S-等位 409 语义一致);**C 近交衰退分析(计算端点不建表)**——`POST /statistics/inbreeding-depression`(入参 trait_id/trait_code + year/trial_study_id 可选 + min_n=10):`blup.relationship_matrix` 的 F=A_ii−1(按树,**非近交树从 individuals 默认 F=0 兜底**否则被排除)+ trait 树级表型 group_by(tree_id) 均值 → **纯 numpy lstsq 回归 phenotype~F**——`sufficient/reason/n_obs/mean_F/max_F/direction/regression{intercept,slope,r2,t_slope,depression_rate,has_depression,interpretation}/high_inbreeding`(F 降序前 10 含 tree_no/value);`has_depression=slope<0 且 t_slope≤−1.645`(单侧 95% 大样本近似);n候选/空候选/λ<0 409;**E MABC 标记辅助回交**——`POST /statistics/mabc-progress`(纯计算不落库):**前景选择**复用 `_mas_hit_map`(五语义、stage 无关、童期幼苗可用)`foreground_pass`=命中≥阈值;**背景基因组恢复率**对 background 面板逐标记候选 vs 轮回亲本「纯合剂量一致」比例(`_pure_homologous`,杂合不计)均值=恢复率 %;**回交代建议**按 `_BC_LADDER`(F1→BC1→…→BC5→自交固定)+ `_MABC_GENERATIONS`:晋级下一代 / 再回交一代 / 恢复无法评估(轮回亲本缺背景基因型)补测后判定 / 前景未通过;输出 per_tree{前景命中/恢复率/bc_generation/reason} + `recommended`(前景通过按恢复率降序)+ summary + warning;门禁空候选/空面板/generation 非法/background_target∉(0,100]/轮回亲本不存在 409。e2e `ocs` 17 断言(Σc/c≥0/λ 梯度 avg_kinship 单调非增 1.0→0.5157→0.3211→0.2685/λ=0 退化=贪心/partial+warning/三 409)+ `mabc` 23 断言(前景 5pass+1fail/恢复复算 100·75·75·0/reason 晋级·再回交·前景未通过/recommended 恢复率降序/童期 C6/轮回亲本缺基因型 None+warning/六 409)+ 回归 molrigor/gwas/mating_kinship/mating_ebv/gblup_cv/statsflow/spatial/busflow/gaps_base/b8/prediction_rigor/audit_fix 全套零回归 + vue-tsc EXIT=0 + HTTP 核验。`_DATA_VERSION` v2.24。详见规格 §8.34 / `桃育种系统统计引擎实施记录.md` §18.12。);㉔ **A–F 六项功能级缺口·批4 分子亲本验证 + 试验设计深度(2026-08-05,用户实证 grep 指出 genotype_call 对 mendelian/parentage 零命中 + trial_design 仅 RCBD,逐项复核确认 + 拍板「平方格子+一般循环回退」,独立 e2e + 全量回归零回归)**——**A 孟德尔/亲本一致性检验**(`GET /bre/genotype_call/mendelian/check`,纯计算不落库):子代树样本基因型 vs 声明亲本(tree.dam_id/sire_id→germplasm)逐标记 union 规则(子代等位集 ⊆ dam∪sire 并集)判兼容,单株兼容率=兼容/信息位点、rate0→has_warning);系谱消费方断开隔离树亲本链(退化为 founder 表型保留:`_build_pedigree`/`run_ablup`/`inbreeding_depression`/`kinship_matrix`/`run_cv`/`_germplasm_pedigree`,GBLUP/ssGBLUP 基因型样本不受影响);**N2 统计引擎全量数值真值校验**:新探针 `e2e_stats_numerics_20260805.py` **ok=101 fail=0** 逐段独立 ground-truth 对照(blup/gxe/block/spatial-aniso/gblup/ssgblup/bivariate/multi/FW/ammi/rcbd/diallel/ocs/UPGMA/ΔG/data-quality/kinship 全验证),探针自身两处 fixture 缺陷修复(低 h² 小样本 r_g 假收敛换高遗传力良态设计、B4.2 硬编码 tree_id),**引擎本体零缺陷**;`_DATA_VERSION` v2.26、`blup.ENGINE_VERSION` 1.5.0。详见规格 §8.37 / `桃育种系统统计引擎实施记录.md` §18.15。);㉖ **批6 规划↔目标↔试验 可追溯性闭环(2026-08-05,用户全代码库 grep 实证 `bre_plan` 彻底孤立根表——零传入外键、`bre_trial` 无 plan_id/target_id、无 project/program 表;拍板「加试验挂目标」最小双链,独立 e2e 探针 17 断言 + 全量回归零回归)**——**三可空 FK 双链闭合**:`bre_target.plan_id` + `bre_trial.plan_id`(目标链 plan→target→cross_combination→tree)+ `bre_trial.target_id`(试验链 plan→trial→trial_study→tree),weld SQL 幂等、target/trial 三件套透传 + get_list 解析 `plan_name`/`target_name`;**父校验**(target/trial create 坏 plan_id/target_id → 409,None 跳过——孤目标/独立试验仍可建)+ **删除守卫**(plan.delete 有 target/trial 子 409、target.delete 有 combo/trial 子 409);**计划可追溯端点** `GET /bre/plan/trace/{id}`(plan 基础 + targets[] 组合株数 + trials[] 试验点株数 + totals,软删过滤、孤目标不入链、不存在 409);前端 target/trial 「所属计划/目标」select + plan 页「追溯」弹窗。验证:`e2e_plan_trace_20260805.py` 17 断言全过 + 全量回归零回归(含审计切面残留检查限定 `entity_type` 根治跨表撞 id 误报)+ vue-tsc EXIT=0 + HTTP 核验;`_DATA_VERSION` v2.27。详见规格 §8.38 / `桃育种系统统计引擎实施记录.md` §18.16。);㉗ **批7 统计引擎全量十模块 golden 复查·系谱键空间错位修复(2026-08-05,用户质疑 GBLUP/ssGBLUP、genetic_corr、type_b、cluster、stability、OCS、MABC、inbreeding_depression、trial_design、mendelian 十新建模块缺逐模块数值校验(golden-test 缺失)、历史上 G–L 静默失效都靠校验暴露,拍板「全量十项逐一复查」)**——复查确认 **8/10 已 golden**(GBLUP/ssGBLUP、cluster、stability、OCS、MABC、inbreeding_depression、trial_design、mendelian 既有探针已具严格数值断言),**type_b 与 genetic_corr 补齐 golden 缺口**——审计深处暴露**系谱键空间错位真 bug**:`run_type_b_heredity`/`correlation_matrix`(mtblup)/`run_genetic_corr` 四个服务路径把 **int `tree_id` 表型键**传给 `mtblup.solve_bivariate`,而 `_build_pedigree` 系谱用 **str `"t{id}"`**,`blup.solve` 把系谱外个体自动并入 base(无关 founder)→ **系谱被静默丢弃**(平衡数据同株跨环境掩盖;**分裂家系数据**无同株跨环境、靠全同胞链桥接暴露——int 键 → r_B/r_g 锁死 ±0.99 平剖面端点假收敛假阴性)。修复:共享 helper `_tree_individual_keys`(int→`"t{id}"` 对齐系谱键空间)收口**四个调用点**(`_genetic_corr_matrix`/`_smith_hazel_index` mtblup 分支/`run_genetic_corr`/`run_type_b_heredity`),`run_ablup` 早已正确不受影响。新探针 `e2e_rg_golden_20260805.py` **ok=21 fail=0**(20 founder×10 全同胞家系×3 同胞=30 株、3 站点、高遗传力 h²≈0.89:P0 诊断 str 键恢复 r_B=0.7+/int 键锁死 −0.99;P1 type_b 分裂家系 r_B(1,2)>0.5 + **n_common=0**(无共有树、系谱桥接生效直接证据);P2 correlation 平衡 r_g(A,B)>0.65/AC,BC<0.4/判别间隙>0.3/h²∈(0,1);P2b 分裂家系 discriminator r_g(S1,S2)>0.5)。验证:`e2e_statsflow`(单家系 fixture 对键 bug 不敏感 15 断言零变化)+ `e2e_advanced_mtblup` 回归通过 + **全量回归零回归**(60 套件,含版本断言同步)+ vue-tsc EXIT=0 + 后端重启 HTTP 核验(type-b-heredity / genetic-corr/run 分裂家系判别经真实端点)+ 备份镜像 `Temp/claude/backup_keyspace_20260805/`;`_DATA_VERSION` v2.28。详见规格 §8.39 / `桃育种系统统计引擎实施记录.md` §18.17;㉘ **回归基建:统一回归 + 高风险模块正式 tc 套件(2026-08-05,用户三点质疑——tc 零覆盖分析模块 / 批7 无 golden / 新模块算偏当正常输出(G/I 静默失效同款隐患)——核实全部成立,拍板「两者都要」;纯测试基建不 bump 版本)**——`run_regression.py`(仓库根)顺序串跑 `Temp/claude/e2e_*.py` + `backend/scripts/test_bre_analysis_*_tc.py` = **63 测试**、每测试独立子进程 log 到 `regression_logs/`、EXIT 汇总、`DPB_PYTHON` 覆盖解释器、utf-8 防 GBK;**4 个正式 tc 套件**(TestClient 真实 API + golden 数值断言,高风险模块必测):`test_bre_analysis_ssgblup_tc.py` ok=14(非基因型株 GBLUP EBV=0 / ssGBLUP EBV≠0——H⁻¹ 拼接判别)、`test_bre_analysis_ocs_tc.py` ok=23(投影梯度收敛/λ 单调/c'Ac/λ=0 退化)、`test_bre_analysis_mabc_tc.py` ok=24(前景/背景恢复率/晋级/童期)、`test_bre_analysis_mendelian_tc.py` ok=9(POST run 落库 + 隔离注册表 + 阈值 0.0 可逆);教训:HTTP schema 校验 → 422 / 服务层 CustomException → 409、tc 须显式 import UserModel、隔离断言时序、启动残留预清理 `_wipe`。全量回归 58/59(gaps_engines 批量卡死单跑 RC=0 确认);备份镜像 `Temp/claude/backup_test_rigor_20260805/`;交付流程:改统计引擎/业务链后端后跑 `python run_regression.py` 全量 EXIT=0。详见规格 §8.40 / `桃育种系统统计引擎实施记录.md` §18.18。);㉙ **回归基建 2:分析模块 tc 补齐(5 套)+ pytest 接入 + 本机 CI stage(2026-08-05,用户对 4 套 tc 认可后拍板补齐 genetic_corr/type_b_heredity/inbreeding_depression/trial_design/run_stability 五个静默失效风险较低模块,v8 验收标准定为「CI 跑通 pytest 全绿」;纯测试基建不 bump 版本,`_DATA_VERSION` 保持 v2.28)**——**补 5 套正式 tc**(`backend/scripts/test_bre_analysis_{genetic_corr,type_b,inbreeding,trial_design,stability}_tc.py`,TestClient 真实 API + `-X utf8` + `_wipe()` 残留预清理):genetic_corr ok=19(平衡 3 性状 r_g(AB)>0.65/间隙>0.3 + **分裂家系判别** r_g(S1,S2)>0.5 且 n_common=0,独立种子 SEED_S=3)、type_b ok=18(3 站点分裂家系 r_B>0.5 且 n_common=0,复刻 tY-first rng 序→0.990)、inbreeding ok=44(方向感知:衰退 slope<−40/无衰退/上升 slope>30 + deep 系谱 F=0.125 + min_n insufficient)、trial_design ok=25(增广对照 block_no=None 不落库 + α-格子平方 λ≤1/非素数告警/一般贪心 + **block_no=rep*100+block DB 读回** + 四 409)、stability ok=16(FW b≈2/1/0 + flag_stable C1=True + AMMI ASV rank + engine_version 落库);**pytest 接入**——`backend/tests/test_stat_analysis_tc.py` 9 套 @pytest.mark.pg + subprocess 串跑真实脚本、skipif 探测真实 PG16+Redis(**protocol=2**——redis-py 7.x 默认 RESP3 对老 WSL redis-server 报 ResponseError,与 app redis_connect 对齐)、env 显式构建剥 conftest 污染解析 .env.dev 可 TC_CI_* 覆盖;pyproject 注册 pg marker、conftest collect_ignore_glob 排除 scripts;**本机 CI stage `run_ci.py`**(仓库根,step 0 探测 PG+Redis 失败指引 / step 1 `pytest -m pg` / step 2 EXIT=0「CI 全绿」,日志 Temp/claude/ci_logs/;不做 GitHub Actions——全工程无 git 仓库);验证 5 套 EXIT=0(122 断言)+ `pytest -m pg` 9 passed(146s)+ skipif 9 skipped 不红 + run_ci EXIT=0 + run_regression glob 自动纳入;备份镜像 `Temp/claude/backup_analysis_tc2_20260805/`;交付流程:v8 验收 = `python run_ci.py` EXIT=0。详见规格 §8.41 / `桃育种系统统计引擎实施记录.md` §18.19。);㉚ **O2–O4 统计方法学广度(2026-08-05,用户系统核对统计引擎能力后确认三项 P2 方法学缺口——幼年–成年遗传相关 / 非线性基因组预测 / 约束选择指数,均零命中、不影响当前业务正确性属方法库广度,拍板「O2–O4 先落地」+ 两决策:① O4 用经典 Kempthorne–Nordskog =0 闭式解(纯 numpy 不迭代)② 前后端同步扩展;`_DATA_VERSION` v2.28 → v2.29、genomic 1.4.0,无新表)**——**O2 幼年–成年遗传相关**:`_values_by_tree_env` 加 env_dim=stage 分支(`coalesce(obs.stage, TraitModel.stage)` 分组)+ `run_type_b_heredity` env_dim 校验扩 stage——发育阶段当"环境"逐对 `mtblup.solve_bivariate` REML 精确估 r_g(>0.5 表示童期可借力成株、早期间接选择有效),复用 `bre_type_b_result` 不建表,前端 Type-B 下拉加 stage;**O3 非线性基因组预测**:`genomic.py` 1.4.0 加 `solve_rrblup`(岭回归逐标记 u=(Xc'Xc+λI)⁻¹Xc'yc,**与 GBLUP 数学对偶**(同数据 EBV 逐位吻合)、输出 marker effects)+ `solve_bayesb`(BayesB Gibbs 可变选择 δ_j~Bernoulli(1−π)、**确定性 seed 可复现**)+ `kfold_cv_genomic` dosage 分派 rrblup/bayesb(实证 CV);`run_gblup` method 扩 rrblup/bayesb + seed 透传落 PredictionModel(method=RRBLUP/BayesB)+EBV 行、`run_cv` GS 分支 KFCV/RRBLUP·KFCV/BayesB;**O4 约束选择指数**:`_smith_hazel_index` 加 restricted(**Kempthorne-Nordskog Lagrange 闭式投影** b=P⁻¹·G·(I−M)·a、M=C'(C·G·P⁻¹·G·C')⁻¹·C·G·P⁻¹·G → 受限性状 **ΔG=0**;对照旧公式 M=C'(C·P⁻¹·C')⁻¹·C·P⁻¹ 实测 ΔG=4.46 不满足证明修正必要)+ economic_weights(绝对尺度不归一化、可与 restricted 叠加、asc 翻转),method 校验扩 zsum/smith_hazel/restricted、落库 SI_RES;**顺手修复 tc 暴露的既有缺口**:`cv_detail` 折明细缺 `cv_result_id` → 任意 GET /cv/{id} 400(既有 tc 未触达静默失效)——补上可查;验证:e2e 探针 19 断言(对偶逐位 1.42e-13 / BayesB 可复现 / 受限 ΔG=0 机器精度 −8.88e-16 / O2 真实 PG r_g=0.6145)+ 正式 tc 第 10 套 `test_bre_analysis_o234_tc.py` ok=41 fail=0(r_g=0.8995 / bayesb 同 seed 两次 EBV 逐位一致 / rrblup·GBLUP EBV 相关 1.000 / b_受限=0.0 / economic desc>0·asc<0 / 三 409)+ 全量回归零回归 + `run_ci.py` EXIT=0;备份镜像 `Temp/claude/backup_o234_20260805/`。详见规格 §8.42 / `桃育种系统统计引擎实施记录.md` §18.20。;㉛ **O1 花期数据消费补强(2026-08-06,用户系统核对花期字段零消费——`bre_germplasm.flowering_period` 等从未被任何统计/决策逻辑消费,核实确认 + 定性 P1 数据消费补强 + 辅助信息(非 S-allele 级硬错误防呆),拍板三档位落地;`_DATA_VERSION` v2.29 → v2.30,无新表)**——① **花期数据建模升级(结构化区间,前置)**:`bre_germplasm` 加双 Date 列 `bloom_start_date`/`bloom_end_date`(weld SQL 幂等 + COMMENT),`flowering_period` 5 档字典保留作粗粒度 fallback;model/schema/service(Excel 导入导出/模板 header 同步)+ 前端 germplasm 表单改两个独立 date 字段直绑(useCrudForm 无提交前 transform → 弃 daterange,与 acquisition_date 同模式零转换)——花期字段首次被消费;② **mating_recommend 花期重叠软警示(不硬阻断)**:`_bloom_overlap`(结构化区间月日年际比较 `_md()` 取 tm_yday 优先;缺区间 fallback 5 档 `_BLOOM_TIER_ORDER` 相邻/同档视为重叠)四态——overlap(重叠可同期授粉,**不进 flags[] 不阻断**)/ store(父本先开花需采粉贮藏)/ offset(母本先开花注意柱头可授期)/ unknown(缺花期数据提示补录),pairs 加 bloom_status/bloom_msg、非重叠三态追加 flags[](前端红色 tag);③ **授粉窗口规划消费 bre_pollen**:新端点 `GET /bre/pollen/window-plan?window_start=&window_end=`(`window_plan` 服务方法,窗口反向 409)——花期种质 `_project_bloom` 结构化区间按年际月日投影到窗口起始年(3-4 月单年窗口假设)、缺区间 fallback `_BLOOM_TIER_SPAN` 5 档代表区间并标注「粗粒度」;**窗口覆盖三态** covered(整段有批次 collect≤bs and expiry≥be)/partial/none(n_need_pollen 计 none);**库存五态** full/expiring(窗口内到期)/entering(窗口内采集)/expired(窗口前失效)/not_started;前端花粉页「授粉窗口规划」卡(daterange + 生成 + 缺粉 alert + 花期覆盖表 + 库存状态表,默认 03-01~04-30)。验证:e2e 探针 `e2e_o1_bloom_20260805.py` 4 段全断言(花期建模双 Date 落库 / 6 种质 15 对四态(区间三类 + 档位 store + 混合对 unknown)/ 窗口覆盖三态 + 档位投影 + 缺粉集合 / 库存五态 / 反向 409)+ vue-tsc EXIT=0 + 全量回归零回归 + `run_ci.py` EXIT=0 + HTTP 核验;备份镜像 `Temp/claude/backup_o1_bloom_20260805/`。详见规格 §8.43 / `桃育种系统统计引擎实施记录.md` §18.22。 +> **状态约定**:✅ 已落地(基础字段) | ⏳ 待建(规划已定,待生成器落地) | ⏳ 待灌(数据) | — 引擎/端点(不建表)。**2026-08-04 追加**:① 统计引擎两轮 P0 已代码落地(G×E + 性状方向标注);② 花粉档案 §3.19 落地(bre_pollen + 授粉窗口校验);③ 选择指数升级为**无性系级聚合**(决策正确性,可靠加权聚合 + 自株退化 + 系级决选);④ 配合力 GCA/SCA 支持**交配设计 `design_type`**(统计严谨:完全双列/部分双列/line×tester(NCII)/NCIII,`combining.solve` 按设计分支双因素模型 + 运行按设计过滤组合);⑤ **统计严谨·MLOps 复现性**(① Smith-Hazel 遗传相关改 **Calo 可靠性校正** r_g=r_EBV/√(rel_i·rel_j);② **PA 落库** pa=√reliability + 批次 accuracy 真 PA;③ **MLOps 溯源** data_version/input_hash/engine_version(同数据重跑哈希一致、改观测必变);④ **germplasm_id 填充**(亲本级 EBV 查询可用);⑤ **砧木字典扩列** 矮化类/砧穗亲和性;⑥ **空间竞争协变量** covariate=competition(Chebyshev 邻域株数,无坐标显式报错));⑥ **组合得失漏斗 v_combination_funnel**(按组合汇总 花→果→种→苗→定植→树→入选 六级派生指标——结实率/出苗率/选择强度滚到组合层,配合力第一手证据,只读端点可查);⑦ **分子育种「建数据能力」**(§3.7 四表落地:bre_marker / bre_genotyping_dataset / bre_genotype_sample(weld 补 sample_name)/ bre_genotype_call,全量 CRUD + Excel 导入导出 + **VCF 导入**(marker/sample 自动建号、GT 编码、幂等)+ **DNA 指纹防错**(无性系混杂检测,同源不一致/异源一致两类标记));⑧ **精修项四连发**(① **DUS/品种保护「数据能力(三表)」** bre_dus_descriptor/test/observation + UPOV TG/53 桃描述符种子 15 条 + 菜单 900230 真实页与 G 段按钮 900800-900808;② **stage 感知** bre_trait.stage(童期/成株)驱动 BLUP/选择指数/决策预览的过滤与警示 + 批次 stage 落库;③ **k-fold 交叉验证外部验证** ABLUP+GXEBLUP 双支持(sire 家系分层留出、逐折容错,bre_cv_result/fold 落库,MLOps 溯源复用);④ **按 site EBV 公平性报告** site 系统性偏差 deviation/flagged + 排名一致性 Spearman + G×E 模式);⑨ **六项能力完善**(① **MT-BLUP 遗传相关**——成对双性状 REML(`mtblup.solve_bivariate`:y 堆叠 Var(u)=G0⊗A,Va/Ve 取单性状 REML、仅对 ρ 黄金分割求极大精确 REML),`run_genetic_corr` 逐对 bivariate 组装 G0 + Higham 半正定投影,落 `bre_genetic_corr_result`(job_type=GENCORR),Smith-Hazel 增 `g_method="mtblup"`(单对不收敛回退 Calo+warning、G 非正定特征值截断、`g_source` 落库可审计);② **GBLUP/ssGBLUP 引擎**——`genomic.py`:VanRaden G method1/2 + MAF 过滤 + blend 岭;`solve_gblup` 仅基因型株、`solve_ssgblup` 单步法 H⁻¹=A⁻¹+[[0,0],[0,G⁻¹−A22⁻¹]] 且 **V 块用 H 协方差**;`run_gblup` 样本 tree 直连/名称兜底、未映射+多等位跳过+warning,落 PredictionModel(GBLUP/ssGBLUP)+EBV 行;③ **DUS 特异性统计检验**——QN LSD(t_crit 由 `fdist.f_icdf` 二分反解 F(1,df))+ PQ/QL 状态众数、**必测描述符**驱动建议、参照=同 trial_study 自动 + 显式覆盖、落 `analysis_json` + conclusion 自动建议;④ **AMMI/Finlay-Wilkinson 稳定性**——FW 环境指数回归 b/R²/se_b/flag_stable(完美拟合 b≈1 判稳)、AMMI 残差 SVD→IPC/ASV/ecovalence、`run_stability` site/year 两维、格子 <2×2→409,落 `bre_stability_result`;⑤ **MLOps 重训+版本回滚**——`bre_prediction.is_active` 版本链(无 active 时首批 active)、`model/drift` 输入哈希漂移检测(轻量 `_gather_digest_inputs`)、`model/retrain` 漂移才重训+active 转移、`model/activate` 回滚原语;⑥ **BLUP stage 拆分**——`bre_trait_observation.stage` 观测级发育阶段、`run_ablup(stage)` 按 coalesce(obs.stage,trait.stage) 过滤+model_name `_{stage}` 后缀+批次 stage 落库、未指定但有分歧→note 警示);⑩ **田间试验精度补强·砧木×接穗随机互作(G×R)**——rootstock 由固定哑变量升级为**第三随机效应槽**(`run_ablup(gxr=True)` → method=GXRBLUP/job_type=GXR,按 `(基因型,砧木)` 分组建 Z₂ 随机互作、克隆坍缩系谱与 G×E 共享、rootstock 强制移出固定效应防共线、不可辨识门禁(同基因型无跨砧木/单元内无重复/残差 df<1)→409、σ²gxr/gxr_ratio/n_cross_rootstock 落批次 note;**Z₂ 单槽位 → G×E/G×R 互斥**(双开 409 + 前端互斥 watcher);求解器 `_solve_gxe` 加 factor_label/factor_name(G×E 警告文案逐字节不变、G×R 显「砧木」)、ENGINE_VERSION 1.2.0;G×E tab 加 G×R 复选框、批次表同显 σ²gxr);⑪ **基因组选择(GS)实证严谨性·GBLUP/ssGBLUP 专用 k-fold 交叉验证**——`run_cv` 加 `dataset_id/method/maf_min`(method=KFCV/GBLUP·KFCV/ssGBLUP),`genomic.kfold_cv_genomic` 掩蔽留出 GEBV vs 表型实证 pearson/RMSE(G/H/Hinv 只建一次、固定种子随机分层——GS 同世代样本无家系树,区别于 ABLUP 的 sire 家系留出)、`_gather_gblup_inputs` 提取(run_gblup↔GS CV 复用)、`_build_h` 缺失基因型补 base、ENGINE_VERSION 1.1.0、_DATA_VERSION v2.11;SSR 多等位 dummy 编码留待后续);⑫ **现代桃育种领域覆盖·S-等位基因交配兼容性 + 抗病/需冷量性状字典**(① `bre_germplasm.s_alleles`(String32,Sf=自交亲和型)+ `cross_combination` create/update 兼容校验——任亲本含 Sf→放行、共享 2 个 S 等位→409「配了不结」硬门禁、共享 1 个→`out.s_compat` 半兼容警示;②③ `bre_trait` 种子 5 条 numeric(细菌性穿孔病/褐腐病/白粉病 0-5 级病情指数 + 需冷量 h + 需热量 GDD,category「抗病性」「生态适应性」,into_ebv=1、direction=asc 进选种目标——统计引擎仅消费 numeric);④ **MAS/QTL/GWAS 闭环**(GWAS 引擎 GLM+PC:SVD 群体结构校正 + 逐标记单标记回归 + Bonferroni/BH-FDR + QTL 区间合并,纯 numpy 复用 fdist 零新依赖;bre_gwas_result/bre_gwas_snp/bre_qtl/bre_mas_panel/bre_mas_panel_marker 五新表;服务端 gwas/run|list|{id} + qtl CRUD + mas-panel CRUD(选标记 favorable_dose/direction/effect);decision_preview 加 marker 条件——童期幼苗无表型/无 EBV 也可被标记辅助选入;前端 900214 GWAS/QTL 页翻转(Manhattan/QQ 裸 echarts + 运行对话框 + 三 tab);e2e 7 段全过 + GBLUP/ssGBLUP 回归零回归,2026-08-04);⑭ **九缺口补齐(2026-08-04,用户重扫 12 项缺口剔除 A1/A4 余 9 项拍板全做,三批独立 e2e + 回归)**——批1 经典侧计算端点(纯计算不建表):**数据质量** `POST /statistics/data-quality`(缺失率/CV + IQR + MAD 稳健 z-score 双法异常标记)+ **遗传增益** `POST /statistics/genetic-gain`(ΔG=k·r_g·σ_A,k 经 Acklam 有理近似 Φ⁻¹ 纯 numpy 零 scipy)+ **主动选配** `POST /statistics/mating-recommend`(S-等位硬过滤 + 亲缘 A 惩罚 score=w_ebv·mid_parent_EBV−w_kin·max(0,r−threshold),不落库);批2 分子侧严谨性:**EMMAX 混合模型 GWAS**(`gwas._emmax`:零模型方差分量复用 `genomic._profile_solve` + eigh 对角化 K→V⁻¹ + 逐标记 GLS t²→fdist,method=emmax,run_gwas 后处理抽共享函数 GLM 逐字节不变)+ **SSR 多等位虚拟编码**(VCF 按 ALT 数定 ssr/snp、`_allelic_from_gt` 每等位一列、`build_g_matrix` 多等位 VanRaden 推广)+ **QTL×E 环境互作**(`POST /statistics/gwas-qtl-x-e` 分环境 GWAS 合并 stable/env-specific/异号 G×QTL 警示,不建表)+ **MAS 面板语义**(bre_mas_panel_marker.mode/favorable_allele/haplotype_group 三列 + `_mas_hit_map` additive/dominance/recessive/allele/haplotype 五语义);批3 求解器核心:**AR1×AR1 空间协方差**(`blup.solve_spatial`:e~N(0,σ²e·R)、R=AR1(ρ)⊗AR1(ρ) R⁻¹ 纯 numpy eigh、REML 对 (h²,ρ) 坐标上升黄金搜索、run_ablup(spatial=True)→method=AR1×AR1、与 G×E/G×R 互斥、缺坐标 409)+ **稀疏 A⁻¹+共轭梯度**(`_build_ainv_sparse` COO + `_cg_solve`、n>N_SPARSE(1000) 阈值分派 solver=sparse-cg、Hutchinson k=100 估 PEV 对角 reliability 近似、稠密路径逐字节不变)。`_DATA_VERSION` v2.13、blup 1.3.0 / genomic 1.2.0 / gwas 1.1.0。详见 `桃育种系统统计引擎实施记录.md`);⑮ **十二项真缺口补齐(2026-08-05,用户全量现状盘点收敛 12 件真缺口、剔除过时「待建」标记,拍板只补这 12 件,分三批独立 e2e + 全量回归零回归)**——批1 引擎三件(solver 扩展独立路径):**ssGWAS**(`gwas._ssgwas`:`build_g_matrix` 岭保正定构 G → `solve_gblup` 全样本 GEBV û → **Wang et al. 单步 GWAS MEM 对角近似**反推标记效应 `ê_j=(M_j'Z'G⁻¹û)/(2Σp_j(1−p_j))`、Wald t=ê/SE → fdist,共享 `_finalize`/`_bh_qvalues`/`_cluster` 后处理,method=ssgwas 透传+前端下拉,gwas 1.2.0)+ **全 MT-BLUP solve_multi**(m 性状 y 堆叠 Var(u)=G0⊗A、**EM-REML 迭代** G0+ve 单性状初始+特征值截断/Higham 投影保正定、max_iter 兜底,输出 G0 全元素+r_g 矩阵+n_iter/converged/warning,`run_genetic_corr(full_mtblup=True)` 一次 m 性状给完整 G0,逐对 bivariate 保活,mtblup 1.1.0)+ **AR1×AR1 双参 ρ aniso**(R_ij=ρ_row^|Δrow|·ρ_col^|Δcol|、ρ_row/ρ_col 交替坐标上升,`run_ablup(spatial_aniso=True)` → method=AR1×AR1(aniso)、note 落 rho_row/rho_col,**aniso=False 与 v1 逐字节一致**,blup 1.4.0);批2 模块四件(标准 CRUD 五件套+前端+菜单翻转+文案清理):**field_operation 农事操作 900051**(op_type 字典施肥/喷药/修剪/灌溉/疏果/套袋/采收 bre_dict 灌入)+ **observation 通用观测 900050**(EAV plot 级通用录入、obs_type numeric/text/date)+ **breeding_report 育种报告 900075**(CRUD + `POST /report/generate` 聚合 bre_prediction/bre_stability_result/bre_cv_result 结构化报告)+ **analysis_dataset 分析数据集 900070**(新表 bre_analysis_dataset jsonb×3 + status,「查看性状值」跳 statistics trait-values 复用既有视图);批3 底座五件:**审计切面**(base_crud 三写注入 CREATE/UPDATE 字段级 diff/DELETE + 批量导入仅汇总 IMPORT 一条 + `GET /bre/audit/list`)+ **统一编号服务**(`app/utils/number_gen.py` NumberGenService + `bre_sequence` + **pg_advisory_xact_lock 原子取号**,迁移 cross_combination/tree/selection_result/seedling 4 处格式保持)+ **超期预警定时任务**(每日 job 按 selection_rule 扫描树 N 年未决策/花粉批次过期/seed_lot 超期 → `bre_alert` + `GET /bre/alert/list` + 前端「预警中心」)+ **备份容灾**(`backend/scripts/pg_backup.sh`/`pg_restore.sh`,Temp/backup 按日期 + N 天轮转,每日 02:00 系统 job)+ **BrAPI + MIAPPE 互操作**(`/brapi/v2` 只读五端点 germplasm/observationvariables/studies/observations/miappe,BrAPI 2.x metadata+pagination+result 包裹,bre_germplasm→germplasmDbId/genus=Prunus/species=persica/subtaxa 映射);`_DATA_VERSION` v2.14、gwas 1.2.0 / mtblup 1.1.0 / blup 1.4.0。详见规格 §8.24–8.26 / `桃育种系统统计引擎实施记录.md` §十七);⑯ **桃育种业务流全面完善(2026-08-05,用户要求「对现有功能,结合桃育种实际业务流程与需求,进行全面测试分析,并完善」,三路审计收敛 10 项真缺口,拍板全部实施 + 两裁决:**导入按钮维持隐藏**(23+ 页显式隐藏非缺口)/ **采收入口走农事操作**(op_type=harvest 结构化产量表单,非独立模块),分三批独立 e2e + 全量回归零回归)**——批A 业务数据链(A1-A6):**果实采收结构化入口**(bre_field_operation 加 yield_kg/fruit_count/avg_fruit_weight/marketable_rate,op_type=harvest 树级采收同步写 bre_trait_observation 4 行、plot 级只留操作记录)+ **产量构成/物候/品质/生长量/描述性状种子 19 条**(bre_yield_phenology_traits.sql:单株产量/果数/均果重/好果率/裂果率 + 盛花期/果实成熟期/落叶期 + 可溶性固形物/可滴定酸/固酸比 + 主干周长/树高/冠幅 + 果形/茸毛/离核/果皮底色/着色类型)+ **种质级观测**(bre_observation.germplasm_id 描述数据入口)+ **授粉→收种→种子处理补链**(seed_lot.pollination_id + seed_treatment.seed_lot_id)+ **EAV 观测校验**(numeric/date parse + valid_min/max 越界 + 同(tree/plot/germplasm,trait,year,date)重复 409)+ **育种阶段自动设置 + 成株性状门禁**(tree stage 按 breeding_generation 推断、童期树录成株性状 409);批B 统计算法(B1-B4):**Type-B 多环境遗传力**(run_type_b_heredity 复用 mtblup.solve_bivariate 把环境当"性状"、env_dim=site/year、落 bre_type_b_result + /type-b-heredity + /type-b/{id})+ **遗传相关 genetic 模式**(correlation mode=genetic 逐对 bivariate 组装 G0,pheno 分支逐字节不变、不落库)+ **选择指数自动权重**(auto_weights=True → default_h2 兜底 0.1、强制 use_h2=False、direction 翻转)+ **UPGMA 层次聚类**(run_cluster 纯 numpy,POST /statistics/cluster 不落库);批C 前端(C1-C4):观测录入增强(observation plot/站点/germplasm 过滤 + germplasm 表单、field_operation harvest 结构化表单)+ 统计可视化 ECharts(correlation 热图 + gblup 可靠性热图/趋势折线 + type-B tab + cluster 树状图)+ combination-funnel 前端入口对话框 + stage_phenotype 容器同步(通用观测/农事操作翻真实页);`_DATA_VERSION` v2.15。详见规格 §8.27 / `桃育种系统统计引擎实施记录.md` §十八);⑰ **遗传链完整性两次修复(2026-08-05,用户先后报告「建种质断链」与「近交惩罚静默关闭」,核实确认 + 拍板落地,e2e 验证通过)**——① **世代闭环修复**:`selection_result._promote_tree_to_clone` 单株晋级到 line/regional_trial/released 建种质时**同步落 `bre_pedigree`**(child_code=品种名、dam/sire=tree.dam/sire 回退组合 female/male_parent_id、combination_id、generation、child_code 幂等)→ 晋升种质不再被统计引擎当 founder,跨世代系谱链闭合(判别点 A[T2,T3]=0→0.125、A[gT1,gA]=0.5);② **近交惩罚「静默关闭」修复**:`mating_recommend` 亲缘矩阵旧靠 `bre_pedigree.child_code` 精确字符串匹配建系谱,名称大小写/空格/代号错配即全空 → 所有配对 r=0 → kin_pen=0 → 近交规避失效且无告警;重写 `_germplasm_pedigree` **多源级联**(① 规范化 strip+casefold 名称反查显式系谱为权威 → ② `tree.germplasm_id` 直连 FK 兜底(多对去重后唯一一致才采纳、歧义回退 founder)→ ③ founder;祖先 BFS 展开 + 去重,修复候选既是 child 又是祖先时系谱重复抛『系谱个体重复』崩溃);`mating_recommend` 新增 **kinship_status(ok/partial/none) + kinship_warning + 逐对『亲缘未解析(r 按 0 计)』旗标**,绝不静默 no-op;设计红线 **FK 只兜底不覆盖显式系谱**(A6 回归 r=0.5 保真)。e2e 14 组断言 + 回归 A6/batch/gaps_base 零回归 + HTTP 核验新字段。详见规格 §8.28 / `桃育种系统统计引擎实施记录.md` §18.5/§18.6。);⑱ **选配「亲本 EBV」回退逻辑方法学修复(2026-08-05,用户报告 `mating_recommend` 旧实现对缺直接 EBV 的候选用「其子代树 EBV 均值」代表亲本育种值——① progeny-test 值非亲本自身值;② 双亲都缺时同一子代树等权计入 dam 与 sire 两侧稀释;③ 跨组合混合无区分,会误导选配排序。核实确认 + 拍板落地)**——`mating_recommend` EBV 回退重写:**direct(germplasm_id 直连亲本自身预测值)优先**、有 direct 绝不被子代树值污染;回退仅当无自身值时用、**显式标注「子代测验近似」并降权 0.5 参与打分**(`_EBV_SOURCE_FACTOR={direct:1.0,progeny:0.5,missing:0.0}`,mid=(f_a·ebv_a+f_b·ebv_b)/2,回退不掩盖 direct);**双亲均缺直接值的子代树剔除**(防双侧等权重复计);记 n_progeny/n_combos 辨别跨组合;返回体加 `candidate_ebv` + `ebv_source.coverage` + 逐对旗标「亲本EBV为子代测验近似(降权0.5)」/「亲本无EBV(按0计)」——EBV 来源显式可查绝不静默降级。e2e 5 组断言 + 回归 A6/亲缘/batch/gaps_base 零回归 + HTTP 核验(真实数据 coverage={direct:0,progeny:0,missing:12})。详见规格 §8.29 / `桃育种系统统计引擎实施记录.md` §18.7。);⑲ **数据质量报告「离群值」方向语义错位修复(2026-08-05,用户报告 §8.23 数据质量端点 IQR/MAD 把任何偏离分布的值都标「离群株」,与育种目标方向错位——正向偏离正是要选的精英单株,核实确认 + 拍板落地)**——`data_quality_report` 方向感知分流:按 `bre_trait.direction`(desc/None=高值优、asc=低值优)判定 `desirable_high`,IQR/MAD 命中的株按偏离侧**双通道分类**(`is_high=v>med`、`is_elite=is_high==desirable_high`):**同向偏离 → `extreme_candidates` 极值候选(精英)**(class=extreme_candidate)/ **反向偏离 → `outliers` 疑似录入错误**(class=outlier,高值优的异常低值、低值优的异常高值),两表互斥;返回体加 `direction` + `extreme_candidates` + `summary.n_extreme_candidates`、`has_warning/n_outliers` 语义收紧为疑似录入错误计数;前端数据质量 tab 改双表(精英 success / 错误 danger)+ 选育方向指示 + 无误报绿色 info alert。e2e 17 株断言(精英 16.0 进 extreme_candidates/side=high、错误 2.0 进 outliers/side=low、互斥不串表)+ 回归 A6/亲缘/batch/gaps_base 零回归 + HTTP 核验(asc 裂果率 bB_DIR5668 与 desc 果重 dA_DIR5668 方向/双表/计数经真实端点正确返回,真实数据无极端值→双表空无误报)。详见规格 §8.30 / `桃育种系统统计引擎实施记录.md` §18.8。);⑳ **ΔG 投影假设声明 + 工程健壮性打磨(2026-08-05,用户报告两点均非功能 bug,核实确认 + 拍板落地,e2e 验证通过)**——① **ΔG 投影假设显式声明**:`genetic_gain` 旧实现 `next_pheno_mean = pheno_mean ± dg` 隐含「环境均值恒定」,年际气候漂移会使投影偏;属投影演示性质,返回体新增 `note` 假设声明(环境均值恒定、未计入年际气候漂移与栽培措施变化;ΔG=k·r_g·σ_A 仅含遗传增益分量,next_pheno_mean 为投影演示值而非实测预测),前端 gain tab info alert 展示;② **空间模型缺坐标优雅降级**:`run_ablup(spatial=True)` 旧实现任一观测树缺 row/col 即整批 ValueError,现改为**剔除缺坐标株表型(保留系谱 → EBV 仍由系谱预测)、仅全缺坐标才整批 409**,剔除数显式写批次 note「缺坐标剔除N株表型(EBV由系谱预测)」、input_hash 用剔除后子集;③ **numpy 工程清理**:7 处函数体内重复 `import numpy` 收敛为模块级单次 + `data_quality_report` 逐株 flag 循环改 numpy 布尔掩码向量化(输出逐字节不变)。e2e spatial 降级场景(5 株含 1 缺坐标 → method=AR1×AR1/train_n=4/note 含剔除/5 株均有 EBV)+ A6 note 断言 + 回归 spatial/A6/亲缘/batch/gaps_base 五套零回归 + vue-tsc EXIT=0 + HTTP 核验(genetic-gain 返回 note、data-quality 返回 direction)。详见规格 §8.31 / `桃育种系统统计引擎实施记录.md` §18.9。);㉑ **A–F 六项功能级缺口·批1 父本验证 + 近交衰退(2026-08-05,用户全模块 grep 实证 A–F 六项缺口均为「有无/深度」问题、非正确性 bug,逐项复核确认 + 拍板「全部六项分批落地」;批1 独立 e2e + 全量回归零回归)**——**A 父本验证(consistency 不落库)**:① `pollination` 花粉来源一致性——`_validate_pollen_parent`(create/update 统一):花粉批次声明采集父本树(pollen.male_tree_id)时授粉声明父本树须一致否则 409(含「不一致」文案),批次无采集父本(混合/外地)或授粉未声明放行不武断;② `tree` 显式亲本 vs 组合父母本警示——`_resolve_parents` 返回 `(dam,sire,warning)`:显式 dam/sire 偏离组合 female/male_parent 时返回体 `warning`「显式母本/父本与组合母本/父本不一致(以显式为准)」(双侧/单侧聚合),**warning 不阻断**(tree 可独立建,武断 409 误伤),schema 只读字段不落库;③ `cross_combination` 同质自交门禁——`_check_parent_roles` 中 `female_parent_id==male_parent_id` → 409「母本与父本同为「{名}」,自交组合请确认」(create/update 生效亲本统一,与 S-等位 409 语义一致);**C 近交衰退分析(计算端点不建表)**——`POST /statistics/inbreeding-depression`(入参 trait_id/trait_code + year/trial_study_id 可选 + min_n=10):`blup.relationship_matrix` 的 F=A_ii−1(按树,**非近交树从 individuals 默认 F=0 兜底**否则被排除)+ trait 树级表型 group_by(tree_id) 均值 → **纯 numpy lstsq 回归 phenotype~F**——`sufficient/reason/n_obs/mean_F/max_F/direction/regression{intercept,slope,r2,t_slope,depression_rate,has_depression,interpretation}/high_inbreeding`(F 降序前 10 含 tree_no/value);`has_depression=slope<0 且 t_slope≤−1.645`(单侧 95% 大样本近似);n候选/空候选/λ<0 409;**E MABC 标记辅助回交**——`POST /statistics/mabc-progress`(纯计算不落库):**前景选择**复用 `_mas_hit_map`(五语义、stage 无关、童期幼苗可用)`foreground_pass`=命中≥阈值;**背景基因组恢复率**对 background 面板逐标记候选 vs 轮回亲本「纯合剂量一致」比例(`_pure_homologous`,杂合不计)均值=恢复率 %;**回交代建议**按 `_BC_LADDER`(F1→BC1→…→BC5→自交固定)+ `_MABC_GENERATIONS`:晋级下一代 / 再回交一代 / 恢复无法评估(轮回亲本缺背景基因型)补测后判定 / 前景未通过;输出 per_tree{前景命中/恢复率/bc_generation/reason} + `recommended`(前景通过按恢复率降序)+ summary + warning;门禁空候选/空面板/generation 非法/background_target∉(0,100]/轮回亲本不存在 409。e2e `ocs` 17 断言(Σc/c≥0/λ 梯度 avg_kinship 单调非增 1.0→0.5157→0.3211→0.2685/λ=0 退化=贪心/partial+warning/三 409)+ `mabc` 23 断言(前景 5pass+1fail/恢复复算 100·75·75·0/reason 晋级·再回交·前景未通过/recommended 恢复率降序/童期 C6/轮回亲本缺基因型 None+warning/六 409)+ 回归 molrigor/gwas/mating_kinship/mating_ebv/gblup_cv/statsflow/spatial/busflow/gaps_base/b8/prediction_rigor/audit_fix 全套零回归 + vue-tsc EXIT=0 + HTTP 核验。`_DATA_VERSION` v2.24。详见规格 §8.34 / `桃育种系统统计引擎实施记录.md` §18.12。);㉔ **A–F 六项功能级缺口·批4 分子亲本验证 + 试验设计深度(2026-08-05,用户实证 grep 指出 genotype_call 对 mendelian/parentage 零命中 + trial_design 仅 RCBD,逐项复核确认 + 拍板「平方格子+一般循环回退」,独立 e2e + 全量回归零回归)**——**A 孟德尔/亲本一致性检验**(`GET /bre/genotype_call/mendelian/check`,纯计算不落库):子代树样本基因型 vs 声明亲本(tree.dam_id/sire_id→germplasm)逐标记 union 规则(子代等位集 ⊆ dam∪sire 并集)判兼容,单株兼容率=兼容/信息位点、rate0→has_warning);系谱消费方断开隔离树亲本链(退化为 founder 表型保留:`_build_pedigree`/`run_ablup`/`inbreeding_depression`/`kinship_matrix`/`run_cv`/`_germplasm_pedigree`,GBLUP/ssGBLUP 基因型样本不受影响);**N2 统计引擎全量数值真值校验**:新探针 `e2e_stats_numerics_20260805.py` **ok=101 fail=0** 逐段独立 ground-truth 对照(blup/gxe/block/spatial-aniso/gblup/ssgblup/bivariate/multi/FW/ammi/rcbd/diallel/ocs/UPGMA/ΔG/data-quality/kinship 全验证),探针自身两处 fixture 缺陷修复(低 h² 小样本 r_g 假收敛换高遗传力良态设计、B4.2 硬编码 tree_id),**引擎本体零缺陷**;`_DATA_VERSION` v2.26、`blup.ENGINE_VERSION` 1.5.0。详见规格 §8.37 / `桃育种系统统计引擎实施记录.md` §18.15。);㉖ **批6 规划↔目标↔试验 可追溯性闭环(2026-08-05,用户全代码库 grep 实证 `bre_plan` 彻底孤立根表——零传入外键、`bre_trial` 无 plan_id/target_id、无 project/program 表;拍板「加试验挂目标」最小双链,独立 e2e 探针 17 断言 + 全量回归零回归)**——**三可空 FK 双链闭合**:`bre_target.plan_id` + `bre_trial.plan_id`(目标链 plan→target→cross_combination→tree)+ `bre_trial.target_id`(试验链 plan→trial→trial_study→tree),weld SQL 幂等、target/trial 三件套透传 + get_list 解析 `plan_name`/`target_name`;**父校验**(target/trial create 坏 plan_id/target_id → 409,None 跳过——孤目标/独立试验仍可建)+ **删除守卫**(plan.delete 有 target/trial 子 409、target.delete 有 combo/trial 子 409);**计划可追溯端点** `GET /bre/plan/trace/{id}`(plan 基础 + targets[] 组合株数 + trials[] 试验点株数 + totals,软删过滤、孤目标不入链、不存在 409);前端 target/trial 「所属计划/目标」select + plan 页「追溯」弹窗。验证:`e2e_plan_trace_20260805.py` 17 断言全过 + 全量回归零回归(含审计切面残留检查限定 `entity_type` 根治跨表撞 id 误报)+ vue-tsc EXIT=0 + HTTP 核验;`_DATA_VERSION` v2.27。详见规格 §8.38 / `桃育种系统统计引擎实施记录.md` §18.16。);㉗ **批7 统计引擎全量十模块 golden 复查·系谱键空间错位修复(2026-08-05,用户质疑 GBLUP/ssGBLUP、genetic_corr、type_b、cluster、stability、OCS、MABC、inbreeding_depression、trial_design、mendelian 十新建模块缺逐模块数值校验(golden-test 缺失)、历史上 G–L 静默失效都靠校验暴露,拍板「全量十项逐一复查」)**——复查确认 **8/10 已 golden**(GBLUP/ssGBLUP、cluster、stability、OCS、MABC、inbreeding_depression、trial_design、mendelian 既有探针已具严格数值断言),**type_b 与 genetic_corr 补齐 golden 缺口**——审计深处暴露**系谱键空间错位真 bug**:`run_type_b_heredity`/`correlation_matrix`(mtblup)/`run_genetic_corr` 四个服务路径把 **int `tree_id` 表型键**传给 `mtblup.solve_bivariate`,而 `_build_pedigree` 系谱用 **str `"t{id}"`**,`blup.solve` 把系谱外个体自动并入 base(无关 founder)→ **系谱被静默丢弃**(平衡数据同株跨环境掩盖;**分裂家系数据**无同株跨环境、靠全同胞链桥接暴露——int 键 → r_B/r_g 锁死 ±0.99 平剖面端点假收敛假阴性)。修复:共享 helper `_tree_individual_keys`(int→`"t{id}"` 对齐系谱键空间)收口**四个调用点**(`_genetic_corr_matrix`/`_smith_hazel_index` mtblup 分支/`run_genetic_corr`/`run_type_b_heredity`),`run_ablup` 早已正确不受影响。新探针 `e2e_rg_golden_20260805.py` **ok=21 fail=0**(20 founder×10 全同胞家系×3 同胞=30 株、3 站点、高遗传力 h²≈0.89:P0 诊断 str 键恢复 r_B=0.7+/int 键锁死 −0.99;P1 type_b 分裂家系 r_B(1,2)>0.5 + **n_common=0**(无共有树、系谱桥接生效直接证据);P2 correlation 平衡 r_g(A,B)>0.65/AC,BC<0.4/判别间隙>0.3/h²∈(0,1);P2b 分裂家系 discriminator r_g(S1,S2)>0.5)。验证:`e2e_statsflow`(单家系 fixture 对键 bug 不敏感 15 断言零变化)+ `e2e_advanced_mtblup` 回归通过 + **全量回归零回归**(60 套件,含版本断言同步)+ vue-tsc EXIT=0 + 后端重启 HTTP 核验(type-b-heredity / genetic-corr/run 分裂家系判别经真实端点)+ 备份镜像 `Temp/claude/backup_keyspace_20260805/`;`_DATA_VERSION` v2.28。详见规格 §8.39 / `桃育种系统统计引擎实施记录.md` §18.17;㉘ **回归基建:统一回归 + 高风险模块正式 tc 套件(2026-08-05,用户三点质疑——tc 零覆盖分析模块 / 批7 无 golden / 新模块算偏当正常输出(G/I 静默失效同款隐患)——核实全部成立,拍板「两者都要」;纯测试基建不 bump 版本)**——`run_regression.py`(仓库根)顺序串跑 `Temp/claude/e2e_*.py` + `backend/scripts/test_bre_analysis_*_tc.py` = **63 测试**、每测试独立子进程 log 到 `regression_logs/`、EXIT 汇总、`DPB_PYTHON` 覆盖解释器、utf-8 防 GBK;**4 个正式 tc 套件**(TestClient 真实 API + golden 数值断言,高风险模块必测):`test_bre_analysis_ssgblup_tc.py` ok=14(非基因型株 GBLUP EBV=0 / ssGBLUP EBV≠0——H⁻¹ 拼接判别)、`test_bre_analysis_ocs_tc.py` ok=23(投影梯度收敛/λ 单调/c'Ac/λ=0 退化)、`test_bre_analysis_mabc_tc.py` ok=24(前景/背景恢复率/晋级/童期)、`test_bre_analysis_mendelian_tc.py` ok=9(POST run 落库 + 隔离注册表 + 阈值 0.0 可逆);教训:HTTP schema 校验 → 422 / 服务层 CustomException → 409、tc 须显式 import UserModel、隔离断言时序、启动残留预清理 `_wipe`。全量回归 58/59(gaps_engines 批量卡死单跑 RC=0 确认);备份镜像 `Temp/claude/backup_test_rigor_20260805/`;交付流程:改统计引擎/业务链后端后跑 `python run_regression.py` 全量 EXIT=0。详见规格 §8.40 / `桃育种系统统计引擎实施记录.md` §18.18。);㉙ **回归基建 2:分析模块 tc 补齐(5 套)+ pytest 接入 + 本机 CI stage(2026-08-05,用户对 4 套 tc 认可后拍板补齐 genetic_corr/type_b_heredity/inbreeding_depression/trial_design/run_stability 五个静默失效风险较低模块,v8 验收标准定为「CI 跑通 pytest 全绿」;纯测试基建不 bump 版本,`_DATA_VERSION` 保持 v2.28)**——**补 5 套正式 tc**(`backend/scripts/test_bre_analysis_{genetic_corr,type_b,inbreeding,trial_design,stability}_tc.py`,TestClient 真实 API + `-X utf8` + `_wipe()` 残留预清理):genetic_corr ok=19(平衡 3 性状 r_g(AB)>0.65/间隙>0.3 + **分裂家系判别** r_g(S1,S2)>0.5 且 n_common=0,独立种子 SEED_S=3)、type_b ok=18(3 站点分裂家系 r_B>0.5 且 n_common=0,复刻 tY-first rng 序→0.990)、inbreeding ok=44(方向感知:衰退 slope<−40/无衰退/上升 slope>30 + deep 系谱 F=0.125 + min_n insufficient)、trial_design ok=25(增广对照 block_no=None 不落库 + α-格子平方 λ≤1/非素数告警/一般贪心 + **block_no=rep*100+block DB 读回** + 四 409)、stability ok=16(FW b≈2/1/0 + flag_stable C1=True + AMMI ASV rank + engine_version 落库);**pytest 接入**——`backend/tests/test_stat_analysis_tc.py` 9 套 @pytest.mark.pg + subprocess 串跑真实脚本、skipif 探测真实 PG16+Redis(**protocol=2**——redis-py 7.x 默认 RESP3 对老 WSL redis-server 报 ResponseError,与 app redis_connect 对齐)、env 显式构建剥 conftest 污染解析 .env.dev 可 TC_CI_* 覆盖;pyproject 注册 pg marker、conftest collect_ignore_glob 排除 scripts;**本机 CI stage `run_ci.py`**(仓库根,step 0 探测 PG+Redis 失败指引 / step 1 `pytest -m pg` / step 2 EXIT=0「CI 全绿」,日志 Temp/claude/ci_logs/;不做 GitHub Actions——全工程无 git 仓库);验证 5 套 EXIT=0(122 断言)+ `pytest -m pg` 9 passed(146s)+ skipif 9 skipped 不红 + run_ci EXIT=0 + run_regression glob 自动纳入;备份镜像 `Temp/claude/backup_analysis_tc2_20260805/`;交付流程:v8 验收 = `python run_ci.py` EXIT=0。详见规格 §8.41 / `桃育种系统统计引擎实施记录.md` §18.19。);㉚ **O2–O4 统计方法学广度(2026-08-05,用户系统核对统计引擎能力后确认三项 P2 方法学缺口——幼年–成年遗传相关 / 非线性基因组预测 / 约束选择指数,均零命中、不影响当前业务正确性属方法库广度,拍板「O2–O4 先落地」+ 两决策:① O4 用经典 Kempthorne–Nordskog =0 闭式解(纯 numpy 不迭代)② 前后端同步扩展;`_DATA_VERSION` v2.28 → v2.29、genomic 1.4.0,无新表)**——**O2 幼年–成年遗传相关**:`_values_by_tree_env` 加 env_dim=stage 分支(`coalesce(obs.stage, TraitModel.stage)` 分组)+ `run_type_b_heredity` env_dim 校验扩 stage——发育阶段当"环境"逐对 `mtblup.solve_bivariate` REML 精确估 r_g(>0.5 表示童期可借力成株、早期间接选择有效),复用 `bre_type_b_result` 不建表,前端 Type-B 下拉加 stage;**O3 非线性基因组预测**:`genomic.py` 1.4.0 加 `solve_rrblup`(岭回归逐标记 u=(Xc'Xc+λI)⁻¹Xc'yc,**与 GBLUP 数学对偶**(同数据 EBV 逐位吻合)、输出 marker effects)+ `solve_bayesb`(BayesB Gibbs 可变选择 δ_j~Bernoulli(1−π)、**确定性 seed 可复现**)+ `kfold_cv_genomic` dosage 分派 rrblup/bayesb(实证 CV);`run_gblup` method 扩 rrblup/bayesb + seed 透传落 PredictionModel(method=RRBLUP/BayesB)+EBV 行、`run_cv` GS 分支 KFCV/RRBLUP·KFCV/BayesB;**O4 约束选择指数**:`_smith_hazel_index` 加 restricted(**Kempthorne-Nordskog Lagrange 闭式投影** b=P⁻¹·G·(I−M)·a、M=C'(C·G·P⁻¹·G·C')⁻¹·C·G·P⁻¹·G → 受限性状 **ΔG=0**;对照旧公式 M=C'(C·P⁻¹·C')⁻¹·C·P⁻¹ 实测 ΔG=4.46 不满足证明修正必要)+ economic_weights(绝对尺度不归一化、可与 restricted 叠加、asc 翻转),method 校验扩 zsum/smith_hazel/restricted、落库 SI_RES;**顺手修复 tc 暴露的既有缺口**:`cv_detail` 折明细缺 `cv_result_id` → 任意 GET /cv/{id} 400(既有 tc 未触达静默失效)——补上可查;验证:e2e 探针 19 断言(对偶逐位 1.42e-13 / BayesB 可复现 / 受限 ΔG=0 机器精度 −8.88e-16 / O2 真实 PG r_g=0.6145)+ 正式 tc 第 10 套 `test_bre_analysis_o234_tc.py` ok=41 fail=0(r_g=0.8995 / bayesb 同 seed 两次 EBV 逐位一致 / rrblup·GBLUP EBV 相关 1.000 / b_受限=0.0 / economic desc>0·asc<0 / 三 409)+ 全量回归零回归 + `run_ci.py` EXIT=0;备份镜像 `Temp/claude/backup_o234_20260805/`。详见规格 §8.42 / `桃育种系统统计引擎实施记录.md` §18.20。;㉛ **O1 花期数据消费补强(2026-08-06,用户系统核对花期字段零消费——`bre_germplasm.flowering_period` 等从未被任何统计/决策逻辑消费,核实确认 + 定性 P1 数据消费补强 + 辅助信息(非 S-allele 级硬错误防呆),拍板三档位落地;`_DATA_VERSION` v2.29 → v2.30,无新表)**——① **花期数据建模升级(结构化区间,前置)**:`bre_germplasm` 加双 Date 列 `bloom_start_date`/`bloom_end_date`(weld SQL 幂等 + COMMENT),`flowering_period` 5 档字典保留作粗粒度 fallback;model/schema/service(Excel 导入导出/模板 header 同步)+ 前端 germplasm 表单改两个独立 date 字段直绑(useCrudForm 无提交前 transform → 弃 daterange,与 acquisition_date 同模式零转换)——花期字段首次被消费;② **mating_recommend 花期重叠软警示(不硬阻断)**:`_bloom_overlap`(结构化区间月日年际比较 `_md()` 取 tm_yday 优先;缺区间 fallback 5 档 `_BLOOM_TIER_ORDER` 相邻/同档视为重叠)四态——overlap(重叠可同期授粉,**不进 flags[] 不阻断**)/ store(父本先开花需采粉贮藏)/ offset(母本先开花注意柱头可授期)/ unknown(缺花期数据提示补录),pairs 加 bloom_status/bloom_msg、非重叠三态追加 flags[](前端红色 tag);③ **授粉窗口规划消费 bre_pollen**:新端点 `GET /bre/pollen/window-plan?window_start=&window_end=`(`window_plan` 服务方法,窗口反向 409)——花期种质 `_project_bloom` 结构化区间按年际月日投影到窗口起始年(3-4 月单年窗口假设)、缺区间 fallback `_BLOOM_TIER_SPAN` 5 档代表区间并标注「粗粒度」;**窗口覆盖三态** covered(整段有批次 collect≤bs and expiry≥be)/partial/none(n_need_pollen 计 none);**库存五态** full/expiring(窗口内到期)/entering(窗口内采集)/expired(窗口前失效)/not_started;前端花粉页「授粉窗口规划」卡(daterange + 生成 + 缺粉 alert + 花期覆盖表 + 库存状态表,默认 03-01~04-30)。验证:e2e 探针 `e2e_o1_bloom_20260805.py` 4 段全断言(花期建模双 Date 落库 / 6 种质 15 对四态(区间三类 + 档位 store + 混合对 unknown)/ 窗口覆盖三态 + 档位投影 + 缺粉集合 / 库存五态 / 反向 409)+ vue-tsc EXIT=0 + 全量回归零回归 + `run_ci.py` EXIT=0 + HTTP 核验;备份镜像 `Temp/claude/backup_o1_bloom_20260805/`。详见规格 §8.43 / `桃育种系统统计引擎实施记录.md` §18.22。;**2026-08-06 追加**:业务链路与数据模型说明书 v1.1 已与代码同步(O1 花期 / N4 血缘 / 亲本校验 / 世代闭环 / §九 QA 与 CI)。 --- diff --git a/doc/桃育种观测业务梳理.docx b/doc/桃育种观测业务梳理.docx new file mode 100644 index 0000000..177f8b0 Binary files /dev/null and b/doc/桃育种观测业务梳理.docx differ diff --git a/frontend/web/src/components/forms/composables/useFormBase.ts b/frontend/web/src/components/forms/composables/useFormBase.ts index 4a69d77..7fd9247 100644 --- a/frontend/web/src/components/forms/composables/useFormBase.ts +++ b/frontend/web/src/components/forms/composables/useFormBase.ts @@ -33,6 +33,9 @@ export interface SanitizeOutputOptions { /** 传递给组件时需排除的表单配置属性 */ const ROOT_PROPS = ["label", "labelWidth", "key", "type", "hidden", "span", "slots"]; +/** 选项超过该数量时自动切换为虚拟滚动选择器(ElSelectV2),避免渲染数万个 el-option 卡死页面 */ +const VIRTUAL_SELECT_THRESHOLD = 2000; + /** 日期选择器类型列表(getProps 中用于传递 type 到 FaDatePicker) */ const DATE_PICKER_TYPES = ["date", "daterange", "datetime", "datetimerange", "monthrange"]; @@ -167,6 +170,18 @@ export const getProps = (item: FormItemBase): Record => { return props; }; +/** + * 大选项集 select 是否切换为虚拟滚动(ElSelectV2)。 + * 有自定义插槽 / render 的表单项保持原样,避免破坏自定义选项渲染。 + */ +export const isVirtualizedSelect = (item: FormItemBase): boolean => { + if (item.type !== "select" || item.slots) return false; + const props = getProps(item); + if (props?.default) return false; + const opts = props?.options; + return Array.isArray(opts) && opts.length > VIRTUAL_SELECT_THRESHOLD; +}; + /** * 获取插槽 —— 过滤掉未定义的插槽 */ diff --git a/frontend/web/src/components/forms/fa-form/index.vue b/frontend/web/src/components/forms/fa-form/index.vue index 48a08f9..722218e 100644 --- a/frontend/web/src/components/forms/fa-form/index.vue +++ b/frontend/web/src/components/forms/fa-form/index.vue @@ -40,8 +40,10 @@ @update:model-value="setFieldValue(item.key, $event)" v-bind="getProps(item)" > - -