周四早会机加车间质量角。这批轴类零件连续做了 12 个批次质量工程师小林把一摞首检报告摊开图纸关键尺寸 Φ40h7公差带 ±0.012mm。前 5 批尺寸稳在 39.992~39.996第 6 批开始往小了飘第 9 批出了两根 39.978直接超下差。客户问是不是设备老化了可我们三台数控同干说不清是哪台、哪批、哪个因素起的头。我点开他们导出的批次测量 CSV。这表里有什么小林问。每条是某批次某件某特征的实际测量值我指着屏幕带 batch_id、设备号、特征名、标称值、实测值。但它只给单件值没按批次算均值/标准差也没做批次间稳定性对比。现在靠翻报告肉眼比12 批翻下来眼都花了还分不出批次自然波动和批次漂移。我就想干一件事小林说把每个批次每个关键尺寸聚合成均值、标准差、CPK再算批次间的波动稳定性哪些批次开始变脸哪些设备做的批次更飘画出来。最好能用控制图思路标出失控批再做个聚类区分稳定批/漂移批/异常批。比如 B06 批均值 39.988标准差 0.0042CPK 1.08B01 批均值 39.994标准差 0.0018CPK 2.1我接话均值往下掉 0.006看着不大但 CPK 从 2.1 掉到 1.08就是过程在退化。再往后 B09 均值 39.981已经出界。对小林点头还有想看是不是某台数控换刀后引起的按设备拆批次看别一股脑算总账。用 pandas 按批次×特征×设备聚合numpy 算均值/标准差/CPK/批次间滑动方差matplotlib 画控制图批次箱线图CPK趋势线热力图scipy 做批次间 ANOVA 与均值漂移检验scikit-learn 聚类批次稳定性画像networkx 建批次-设备-特征关系网我开工程数据自包含合成一批多批次多设备加工数据下载就能跑。敲了行原型batch_stat df.groupby([batch_id,feature,machine])[dev].agg([mean,std])cpk np.minimum(usl-mean, mean-lsl) / (3*std)# 控制图判异: 超出 ±3σ 或连续趋势完整版 OOP 封好我说加载器、批次统计器、过程能力器、稳定性分析器、显著性检验器、批次聚类器、关系网、出图器输出稳定批次清单 5 图 报告存 results/。小林凑近看那以后看报告12 批 × 3 设备 × 关键尺寸B01~B05 稳定绿区B06 起进黄区B09/B10 红区超差控制图标出 B06 均值出警戒线、B09 出控制线聚类分稳定批/缓漂移批/突变批/异常批按设备拆发现是 M02 换刀后批次漂移最明显。工艺单直接挂M02 第6批起复核刀补。对我接话批次稳定性不是看单件合不合格是看过程有没有在悄悄退化。数字孪生里建加工过程质量模型这些批次序列就是标定基准。一、实际应用场景真实痛点场景设定多批次轴类/盘类零件在数控产线加工每批首检抽检留测量数据。质量组需要量化批次间尺寸波动识别过程退化而非单件超差定位到具体设备与批次拐点输出可执行的批次稳定性管控清单。现场原话叙事化不是我们不出合格品小林说是过程在悄悄变脸。前几批 CPK 2.0 以上后几批掉到 1.1单件还合格但已经离悬崖不远了。以前等出超差件才处理现在想提前看到均值在挪、标准差在胖。还有设备的事小林补充三台数控混着做总账看都还行拆开发现 M02 做的批次从第6批开始标准差明显变胖是那次换刀后刀补没回零。不按设备拆这锅又背给设备老化了。核心矛盾单件测量流水 与 批次级统计 过程能力演化 控制图判异 批次间稳定性对比 设备维度下钻 稳定性聚类 之间的断层。二、痛点分析映射到滨州职业学院《先进制造技术》课程模型《先进制造技术》模块 本篇痛点对应数控加工与CAD/CAM技术工序尺寸控制、刀补、加工精度 批次尺寸波动 刀补漂移定位先进制造技术基础公差配合、统计过程控制(SPC)、过程能力 CPK演化 控制图判异FMS与先进生产管理多设备协同、批次追溯、质量管控 按设备拆批次 稳定批清单智能制造与数字孪生加工过程质量数字模型 批次序列作孪生标定数据先进制造新模式数据驱动质量预测 稳定性聚类→预测退化批次一句话总结我们需要一个数控多批次加工尺寸波动与稳定性识别程序用pandas 按批次×特征×设备聚合numpy 算均值/标准差/CPK/滑动方差matplotlib 画控制图/箱线图/CPK趋势/热力图scipy 做 ANOVA 与均值漂移检验scikit-learn 聚类批次稳定性画像networkx 建批次-设备-特征关系网实现从单件测量值到批次稳定性分级 退化拐点定位 设备责任定位。三、核心逻辑讲解大白话3.1 问题本质把批次想成同一道菜连做12天把同程序加工想成食堂师傅按同一菜谱做红烧肉连做12天* 单件实测值 每一碗的咸淡* 批次均值 这一天所有碗的平均咸度* 批次标准差 这一天各碗咸淡散得开不开胖不胖* CPK 这道菜稳不稳得住标准口味* 稳定批 每天均值都贴着标称标准差很瘦* 漂移批 均值悄悄往淡了走单看还行连起来就偏了* 异常批 某天突然出了一碗巨淡的超下差* 控制图 画一条中心线上下警戒线谁出线谁有问题* 按设备拆 看是师傅A还是师傅B的手在变3.2 业务逻辑 → 代码映射导入多批次测量数据│▼ BatchLoader (pandas)读取 CSVbatch_id, machine, feature, nominal, measured, part_id算 dev measured - nominal│▼ BatchStatCalculator (numpy/pandas)按 批次×特征×设备 聚合mean / std / n / min / max│▼ CapabilityTracker (numpy)过程能力演化Cp / Cpk 按批次算│▼ StabilityAnalyzer (numpy)稳定性分析批次均值滑动窗口方差控制图规则(±3σ出界, 连续5点同侧趋势)退化拐点检测│▼ BatchSignificanceTester (scipy)批次间对比多批次 ANOVA前段vs后段 t 检验(漂移判定)设备间批次方差 Levene 检验│▼ BatchProfiler (scikit-learn)批次稳定性聚类特征[均值偏移, 标准差, Cpk, 滑动方差]KMeans → 稳定批/缓漂移批/突变批/异常批│▼ BatchGraph (networkx)批次-设备-特征关系网节点批次/设备/特征边权绝对偏移量, 超界标红│▼ BatchVisualizer (matplotlib)可视化1. 批次控制图(含±3σ控制线/±2σ警戒线)2. 各批次箱线图3. CPK趋势折线图4. 设备×批次 均值偏移热力图5. 批次稳定性聚类散点6. 批次-设备关系网│▼ SyntheticBatchGenerator (numpy)合成数据多批次×多设备×关键特征含前稳后漂、换刀突变、设备偏置3.3 为什么不能只看单件是否合格视角 问题单件在公差内 掩盖均值已漂移 0.006看总均值 前5批拉平后看不出退化批次均值标准差CPK演化 提前2~3批看到退化控制图滑动方差 自动标出拐点批次按设备拆 定位到 M02 换刀责任3.4 分析前后对比维度 翻首检报告 本程序批次均值/标准差 手算 自动聚合过程能力 不计算 按批算CPK趋势退化识别 出超差才知 控制图提前预警设备责任 混算 设备维度下钻批次分级 老师傅经验 聚类出四类画像四、OOP 代码实现4.1 项目结构batch_stability_analysis/├── batch_stability_analysis/│ ├── __init__.py│ ├── batch_loader.py # 数据加载│ ├── batch_stat_calculator.py # 批次统计│ ├── capability_tracker.py # CPK演化│ ├── stability_analyzer.py # 控制图/滑动方差│ ├── batch_significance.py # 统计检验(scipy)│ ├── batch_profiler.py # 稳定性聚类(sklearn)│ ├── batch_graph.py # 批次-设备关系网(networkx)│ ├── visualizer.py # 可视化│ └── synthetic_data.py # 合成数据├── tests/│ ├── __init__.py│ └── test_batch_stability.py├── results/│ ├── control_chart.png│ ├── batch_boxplot.png│ ├── cpk_trend.png│ ├── machine_heatmap.png│ ├── stability_cluster.png│ ├── batch_network.png│ ├── batch_stats.csv│ ├── cpk_series.csv│ ├── stability_clusters.csv│ └── batch_report.txt└── run_batch_analysis.py4.2 核心源码detailssummary/summary数控多批次加工测量数据加载器import pandas as pdfrom pathlib import Pathfrom typing import Optionalclass BatchLoader:加载多批次测量CSVdef __init__(self, filepath: str batch_measure.csv,encoding: str utf-8):self.filepath Path(filepath)self.encoding encodingself._raw: Optional[pd.DataFrame] Nonedef load(self) - pd.DataFrame:if not self.filepath.exists():raise FileNotFoundError(f文件不存在: {self.filepath})self._raw pd.read_csv(self.filepath, encodingself.encoding)rename {}for tgt, al in {batch_id: [batch_id, 批次, batch],machine: [machine, 设备号, cnc],feature: [feature, 特征, dim],nominal: [nominal, 标称值, nom],measured: [measured, 实测值, act],part_id: [part_id, 件号, pid],}.items():if tgt not in self._raw.columns:for a in al:if a in self._raw.columns:rename[a] tgtbreakself._raw self._raw.rename(columnsrename)req [batch_id, feature, nominal, measured]miss [c for c in req if c not in self._raw.columns]if miss:raise ValueError(f缺少必要列: {miss})self._raw[batch_id] self._raw[batch_id].astype(str).str.strip()self._raw[machine] self._raw.get(machine, M-01).astype(str).str.strip()self._raw[feature] self._raw[feature].astype(str).str.strip()self._raw[nominal] pd.to_numeric(self._raw[nominal], errorscoerce)self._raw[measured] pd.to_numeric(self._raw[measured], errorscoerce)self._raw self._raw.dropna(subset[nominal, measured]).copy()self._raw[dev] (self._raw[measured] - self._raw[nominal]).round(4)# 批次序号化, 便于趋势分析bset sorted(self._raw[batch_id].unique(),keylambda x: int(.join(filter(str.isdigit, x)) or 0))self._raw[batch_seq] self._raw[batch_id].map({b: i1 for i, b in enumerate(bset)})return self._raw.reset_index(dropTrue)/detailsdetailssummary/summary批次统计 (numpy/pandas)import numpy as npimport pandas as pdfrom typing import Optionalclass BatchStatCalculator:按 批次×特征×设备 聚合def __init__(self):passdef summarize(self, df: pd.DataFrame) - pd.DataFrame:rows []for (b, feat, m), g in df.groupby([batch_id, feature, machine]):dev g[dev].valuesrows.append({batch_id: b,batch_seq: int(g[batch_seq].iloc[0]),feature: feat,machine: m,n: len(dev),mean_dev: round(float(np.mean(dev)), 4),std_dev: round(float(np.std(dev, ddof1)), 4) if len(dev) 1 else 0.0,abs_mean_dev: round(float(np.mean(np.abs(dev))), 4),max_abs_dev: round(float(np.max(np.abs(dev))), 4),min_dev: round(float(np.min(dev)), 4),max_dev: round(float(np.max(dev)), 4),})out pd.DataFrame(rows)return out.sort_values([feature, batch_seq, machine]).reset_index(dropTrue)def by_batch_feature(self, df: pd.DataFrame) - pd.DataFrame:忽略设备, 按批次×特征聚合(总账视图)rows []for (b, feat), g in df.groupby([batch_id, feature]):dev g[dev].valuesrows.append({batch_id: b,batch_seq: int(g[batch_seq].iloc[0]),feature: feat,n: len(dev),mean_dev: round(float(np.mean(dev)), 4),std_dev: round(float(np.std(dev, ddof1)), 4) if len(dev) 1 else 0.0,max_abs_dev: round(float(np.max(np.abs(dev))), 4),})return pd.DataFrame(rows).sort_values([feature, batch_seq]).reset_index(dropTrue)/detailsdetailssummary/summary过程能力演化 CPK (numpy)import numpy as npimport pandas as pdfrom typing import Optionalclass CapabilityTracker:对称公差 ±tol 下:Cpk min(USL-μ, μ-LSL) / (3σ)按批次计算, 形成演化序列def __init__(self, tol: float 0.012):self.tol toldef track(self, stat_df: pd.DataFrame) - pd.DataFrame:out stat_df.copy()usl, lsl self.tol, -self.tolmu out[mean_dev]sigma out[std_dev].replace(0, 1e-6)cpu (usl - mu) / (3 * sigma)cpl (mu - lsl) / (3 * sigma)cpk np.minimum(cpu, cpl)cp (usl - lsl) / (6 * sigma)out[Cp] cp.round(3)out[Cpk] cpk.round(3)out[cap_level] out[Cpk].apply(self._lvl)return outstaticmethoddef _lvl(cpk: float) - str:if cpk 1.33:return A级if cpk 1.0:return B级if cpk 0.67:return C级return D级/detailsdetailssummary/summary稳定性分析: 控制图 滑动方差 拐点 (numpy)import numpy as npimport pandas as pdfrom typing import Optionalclass StabilityAnalyzer:控制图规则(简化):1. 单批均值超出 ±3σ(全局基准) - 失控2. 连续5批同侧偏移 - 趋势漂移3. 滑动窗口方差突增 - 变胖def __init__(self, warn_k: float 2.0, ctrl_k: float 3.0,trend_n: int 5, window: int 3):self.warn_k warn_kself.ctrl_k ctrl_kself.trend_n trend_nself.window windowdef analyze(self, stat_df: pd.DataFrame) - pd.DataFrame:输入为按批次×特征聚合的表(可含设备维度)out stat_df.copy()# 以全局dev分布为基准建控制限all_mean np.mean(out[mean_dev])all_std np.std(out[mean_dev], ddof1) if len(out) 1 else 0.001out[center] round(all_mean, 4)out[warn_up] round(all_mean self.warn_k * all_std, 4)out[warn_lo] round(all_mean - self.warn_k * all_std, 4)out[ctrl_up] round(all_mean self.ctrl_k * all_std, 4)out[ctrl_lo] round(all_mean - self.ctrl_k * all_std, 4)out out.sort_values(batch_seq)flags []rolling_var []vals out[mean_dev].valuesfor i in range(len(out)):v vals[i]if v out[ctrl_up].iloc[i] or v out[ctrl_lo].iloc[i]:f 失控elif v out[warn_up].iloc[i] or v out[warn_lo].iloc[i]:f 警戒else:f 正常flags.append(f)lo max(0, i - self.window 1)seg vals[lo:i1]rolling_var.append(round(float(np.var(seg, ddof1)), 6) if len(seg) 1 else 0.0)out[ctrl_status] flagsout[rolling_var] rolling_var# 连续同侧趋势out[trend_drift] self._trend_flag(out)return out.reset_index(dropTrue)def _trend_flag(self, df: pd.DataFrame) - list:vals df[mean_dev].valuescenter df[center].iloc[0]res [False] * len(df)for i in range(len(vals)):if i self.trend_n - 1:continueseg vals[i-self.trend_n1:i1]if np.all(seg center) or np.all(seg center):res[i] Truereturn res/detailsdetailssummary/summary批次间显著性检验 (scipy)import numpy as npimport pandas as pdfrom scipy import statsfrom typing import Dict, Listclass BatchSignificanceTester:ANOVA / 前后段t检验 / 设备方差Levenedef __init__(self, alpha: float 0.05):self.alpha alphadef anova_by_feature(self, df: pd.DataFrame) - pd.DataFrame:rows []for feat, g in df.groupby(feature):groups [v[dev].values for _, v in g.groupby(batch_id)]groups [x for x in groups if len(x) 1]if len(groups) 2:continuef, p stats.f_oneway(*groups)rows.append({feature: feat,f_statistic: round(float(f), 4),p_value: round(float(p), 4),significant: bool(p self.alpha),})return pd.DataFrame(rows).sort_values(p_value).reset_index(dropTrue)def front_vs_back(self, df: pd.DataFrame, feature: str,split_seq: int) - Dict:前段vs后段均值漂移检验sub df[df[feature] feature]front sub[sub[batch_seq] split_seq][dev].valuesback sub[sub[batch_seq] split_seq][dev].valuesif len(front) 2 or len(back) 2:return {feature: feature, p_value: np.nan, drift: False}stat, p stats.ttest_ind(front, back, equal_varFalse)return {feature: feature,mean_front: round(float(np.mean(front)), 4),mean_back: round(float(np.mean(back)), 4),p_value: round(float(p), 4),drift: bool(p self.alpha),}def machine_var_levene(self, df: pd.DataFrame, feature: str) - Dict:设备间批次标准差齐性检验sub df[df[feature] feature]groups [v[dev].values for _, v in sub.groupby(machine)]groups [x for x in groups if len(x) 1]if len(groups) 2:return {feature: feature, p_value: np.nan}w, p stats.levene(*groups, centermedian)return {feature: feature, levene_p: round(float(p), 4),var_differs: bool(p self.alpha)}/detailsdetailssummary/summary批次稳定性聚类 (scikit-learn)import numpy as npimport pandas as pdfrom sklearn.cluster import KMeansfrom sklearn.preprocessing import StandardScalerfrom typing import Optionalclass BatchProfiler:基于均值偏移/标准差/CPK/滑动方差聚类def __init__(self, n_clusters: int 4, random_state: int 42):self.n_clusters n_clustersself.random_state random_statedef profile(self, df: pd.DataFrame) - pd.DataFrame:data df.copy()# 用批次×特征总账视图if machine in data.columns:data data.drop(columns[machine])data data.sort_values([feature, batch_seq])feat_cols [abs_mean_dev, std_dev, Cpk, rolling_var]for c in feat_cols:if c not in data.columns:data[c] 0.0X data[feat_cols].fillna(0).values.copy()X[:, 2] -X[:, 2] # Cpk反向scaler StandardScaler()Xs scaler.fit_transform(X)km KMeans(n_clustersself.n_clusters, random_stateself.random_state)data[cluster] km.fit_predict(Xs)centers scaler.inverse_transform(km.cluster_centers_)label_map {}for i, c in enumerate(centers):abs_mean, std, neg_cpk, rvar ccpk -neg_cpkif cpk 1.0:label_map[i] 异常批elif abs_mean 0.004 and std 0.003:label_map[i] 突变批elif abs_mean 0.002 or (std 0.0025):label_map[i] 缓漂移批else:label_map[i] 稳定批data[stability_pattern] data[cluster].map(label_map)return datadef pattern_summary(self, df: pd.DataFrame) - pd.DataFrame:if stability_pattern not in df.columns:return pd.DataFrame()rows []for pat, g in df.groupby(stability_pattern):rows.append({pattern: pat,count: len(g),avg_cpk: round(g[Cpk].mean(), 3),avg_abs_mean_dev: round(g[abs_mean_dev].mean(), 4),batches: , .join(g[batch_id].tolist()[:12]),})order {稳定批:0, 缓漂移批:1, 突变批:2, 异常批:3}return pd.DataFrame(rows).sort_values(pattern, keylambda s: s.map(order)).reset_index(dropTrue)/detailsdetailssummary/summary批次-设备-特征关系网 (networkx)import networkx as nximport pandas as pdfrom typing import Optionalclass BatchGraph:建 批次-设备-特征 关系网def __init__(self):self.G nx.DiGraph()def build(self, stat_df: pd.DataFrame,threshold: float 0.012) - nx.DiGraph:self.G.clear()for b in stat_df[batch_id].unique():self.G.add_node(fB:{b}, ntypebatch)for m in stat_df[machine].unique():self.G.add_node(fM:{m}, ntypemachine)for f in stat_df[feature].unique():self.G.add_node(fF:{f}, ntypefeature)for _, r in stat_df.iterrows():bn, mn, fn fB:{r[batch_id]}, fM:{r[machine]}, fF:{r[feature]}w float(r[abs_mean_dev]) * 100self.G.add_edge(mn, bn, weightround(w, 2))self.G.add_edge(bn, fn, weightround(w, 2))if r.get(max_abs_dev, 0) threshold:self.G.edges[bn, fn][over] Truereturn self.Gdef red_items(self) - pd.DataFrame:rows []for u, v, d in self.G.edges(dataTrue):if d.get(over, False):rows.append({edge: f{u}-{v}, weight: d[weight]})return pd.DataFrame(rows).sort_values(weight, ascendingFalse).reset_index(dropTrue)/detailsdetailssummary/summary可视化 (matplotlib)import numpy as npimport pandas as pdimport matplotlib.pyplot as pltfrom pathlib import Pathplt.rcParams[font.sans-serif] [SimHei, DejaVu Sans]plt.rcParams[axes.unicode_minus] Falseclass BatchVisualizer:def __init__(self, results_dir: str results):self.results_dir Path(results_dir)self.results_dir.mkdir(exist_okTrue)def control_chart(self, df, feature: str):控制图(按批次×特征总账)sub df[df[feature] feature].sort_values(batch_seq)fig, ax plt.subplots(figsize(13, 6))x sub[batch_seq].valuesax.plot(x, sub[mean_dev], -o, color#2C3E50,label批次均值偏差, lw1.5, ms4)ax.axhline(sub[center].iloc[0], colorgreen, ls-, lw1, label中心线)ax.axhline(sub[warn_up].iloc[0], colororange, ls--, lw1)ax.axhline(sub[warn_lo].iloc[0], colororange, ls--, lw1, label±2σ警戒)ax.axhline(sub[ctrl_up].iloc[0], colorred, ls--, lw1.2)ax.axhline(sub[ctrl_lo].iloc[0], colorred, ls--, lw1.2, label±3σ控制线)# 标出异常点alarm sub[sub[ctrl_status].isin([失控, 警戒])]ax.scatter(alarm[batch_seq], alarm[mean_dev],c[red if s失控 else orange for s in alarm[ctrl_status]],s80, zorder5, label判异点)for _, r in alarm.iterrows():ax.annotate(r[batch_id], (r[batch_seq], r[mean_dev]),fontsize8, hacenter, vabottom)ax.set_xlabel(批次序号)ax.set_ylabel(均值偏差 (mm))ax.set_title(f批次控制图 - {feature}, fontsize14, fontweightbold)ax.legend(fontsize8, ncol2)ax.grid(alpha0.3)plt.tight_layout()plt.savefig(self.re利用AI解决实际问题如果你觉得这个工具好用欢迎关注长安牧笛