简介针对短程赛跑中运动员速度变化问题这份数学建模文档基于Keller赛跑模型建立动态优化模型重点分析生理条件、内外阻力与冲力限制对速度的影响通过微分方程推导速度与路程随时间变化的表达式并配合MATLAB非线性拟合求解最大速度时刻。文档共1个doc文件大小4.69MB包含问题重述、模型假设、公式推导、实际数据对比与结论评价等完整学习资源尤其给出了某届奥运会百米决赛前6名速度数据的理论值与实测值比较误差较小适合数学建模备赛者或运动生物力学初学者参考。该资源已有126人学习文档结构清晰图表与数据表格完整能帮助读者系统理解从微分方程建模、曲线拟合到结果检验的完整流程并可直接借鉴Keller模型的建模思路用于类似运动情景分析。1. 短程赛跑的速度曲线远不是三段论那么整齐短程赛跑的速度曲线不是什么「先加速、后匀速、再减速」的三段论形状。用 30 Hz 以上频响的激光测速仪去测 100 m 冲刺数据点会以一种很不规整的方式抖动起跑后的力反冲、途中跑的步频波动、终点前的小幅掉速全部叠加在一条只有 812 m/s 的曲线上。数学建模要做的事不是把曲线平滑得好看而是用带物理约束的参数方程把「加速快慢」「峰值速度」「掉速幅度」这几个可比较的数字提取出来。这类工作直接服务于运动表现分析比如判断某次试跑是起跑慢了还是后程动力不足。下面按一次完整的数据处理流程展开先选模型再讲拟合再做多次试跑的横向对比最后落到排错与脚本复用。适合接触过 Python 和 numpy、想用曲线拟合处理时序数据的工程人员读。2. 短跑速度-时间模型的 3 个候选方程先决定拟合曲线的「骨架」给运动员速度变化建模第一步不是写拟合代码而是确定曲线的函数形式。曲线方程决定了后面所有参数的物理含义也决定了拟合失败时的排查方向。这里从我做这类数据的经验出发比较三个常用候选单指数上升、指数上升加高斯衰减、非参数样条。2.1 单指数上升起点简单但兜不住后程掉速单指数模型形式上是最简洁的v(t) v_max · (1 - e^(-t / τ))这个方程描述的是「速度从零开始以时间常数 τ 逼近极限速度 v_max」的过程。起跑后的加速段表现得很好但真实短跑在 4070 m 之后速度会从峰值缓慢下降单指数模型缺少这一项拟合出的 v_max 会被后程数据拉低而 τ 会被起跑段主导形成一种两边都不讨好的妥协。所以它更适合做基线模型用来和更复杂的模型做信息准则对比而不是直接作为最终报告用的方程。2.2 指数上升 高斯衰减兼顾加速段与掉速段的实用复合形式实际处理 30 Hz 测速数据时我常用的是一个带衰减项的复合模型v(t) v_max · (1 - e^(-t / τ)) · (1 - γ · e^(-((t - t_p) / σ)^2))这里第一个括号是指数上升部分控制起跑加速的节奏第二个括号是一个反向高斯钟形以 t_p 为中心、σ 为宽度乘上幅度系数 γ 描述后程速度衰减。当 γ0 时退化为单指数模型所以两个模型是嵌套的可以用似然比检验判断衰减项是否显著。这个函数形式并不是哪本标准教材里的命名而是工程上用来同时描述「加速」和「掉速」的复合写法胜在参数非常直观v_max 是理论峰值速度τ 越小加速越快γ 是掉速幅度t_p 是峰值速度出现的时间。2.3 PCHIP 样条不预设形态但参数不可解释如果完全不想预设曲线的数学形态可以用 PCHIP分段三次 Hermite 插值做非参数化拟合。PCHIP 比普通三次样条好在不会在数据跳变处产生过冲速度曲线起跑段那种陡峭上升它也能兜住。但问题在于它的输出是一组插值系数而不是几个可读的参数没办法直接说「这个运动员的加速时间常数是 1.8 s」。用它的正确姿势是作为参考曲线先用 PCHIP 找到峰值速度和峰值时间的近似值再把这些值作为参数模型的初值。2.4 三个候选的取舍与对比模型参数数量可解释性后程掉速拟合稳定性主要用途单指数上升2强不支持高基线对比、起跑段单独分析指数上升 高斯衰减5强支持中全程速度曲线主模型PCHIP 样条随节点数无隐含高峰值位置估计、初值计算这三个模型不是竞争关系而是配合关系。实际流程是先用样条估计峰值位置再用单指数做基线最后用复合模型做正式拟合通过 AIC 确认衰减项的加入是否值得。下面这段代码把三个模型放在同一张图上对比各自的拟合效果方便直观感受差异import numpy as np import matplotlib.pyplot as plt from scipy.interpolate import PCHIP from scipy.optimize import curve_fit def model_single(t, v_max, tau): return v_max * (1 - np.exp(-t / tau)) def model_composite(t, v_max, tau, gamma, t_p, sigma): rise 1 - np.exp(-t / tau) decay 1 - gamma * np.exp(-((t - t_p) / sigma) ** 2) return v_max * rise * decay # 模拟一段 10 秒、30 Hz 的速度数据 fs 30.0 t np.linspace(0, 10, 301) v_true model_composite(t, 11.0, 1.8, 0.08, 6.0, 3.0) rng np.random.default_rng(4) v_obs v_true rng.normal(0, 0.15, sizet.shape) # 单指数拟合 p_single, _ curve_fit(model_single, t, v_obs, p0[10, 2]) # 复合模型拟合 p_comp, _ curve_fit(model_composite, t, v_obs, p0[10, 2, 0.05, 6, 3], maxfev20000) # PCHIP 插值不做平滑作为参考 pchip PCHIP(t, v_obs) t_fine np.linspace(0, 10, 601) plt.figure(figsize(8, 4)) plt.plot(t, v_obs, ., markersize2, alpha0.5, labelobserved) plt.plot(t_fine, model_single(t_fine, *p_single), labelsingle exp) plt.plot(t_fine, model_composite(t_fine, *p_comp), labelexp decay) plt.plot(t_fine, pchip(t_fine), labelPCHIP) plt.xlabel(time (s)); plt.ylabel(velocity (m/s)) plt.legend()代码里复合模型初值中 t_p 给 6 s是因为真实数据的峰值通常出现在 58 s 区间γ 初值给到 0.05 表示先假设轻微掉速。PCHIP 不需要初值它的作用是在参数模型拟合失败时提供一条不坏的可视化曲线帮助你判断是数据问题还是模型问题。3. 用 scipy 拟合运动员速度曲线预处理、初值设定与边界约束模型定下来之后拟合本身反而容易出问题因为测速数据不可能像上一节模拟数据那样干净。真实数据里有起跑的反冲尖峰、有步频带来的周期性波动、还有传感器偶尔丢帧。这一章先处理数据再进入拟合最后看拟合质量如何检验。3.1 速度测量数据的预处理降噪、起点对齐与重采样激光测速仪输出的原始速度序列噪声相当大尤其是起跑后前 2 秒运动员脚掌蹬离起跑器时的反冲会在曲线上留下一个明显的毛刺。直接用原始数据拟合模型会把毛刺当做一个真实的加速特征τ 的估计值会被严重拉偏。我一般先用 Savitzky-Golay 滤波做一次平滑窗口长度取奇数大约是 0.30.5 s 对应的采样点数from scipy.signal import savgol_filter fs 30.0 # 采样率单位 Hz window int(0.4 * fs) | 1 # 保证奇数0.4s 窗口 v_smooth savgol_filter(v_obs, window_lengthwindow, polyorder3) # 判定起跑点速度第一次持续超过 0.3 m/s 的时刻 mask v_smooth 0.3 idx_start np.argmax(mask) t t - t[idx_start] # 时间轴从起跑点重新归零 v v_smooth[idx_start:]这里把时间零点放在起跑瞬间很重要因为复合模型的时间 t 是「起跑后经过的时间」如果时间轴从仪器开机开始算τ 和 t_p 两个参数的物理含义都会错位。0.3 m/s 这个阈值不是硬性标准它的作用是过滤掉传感器零点漂移导致的虚假起始段如果运动员起跑前有准备姿势的微小移动可以把阈值提高到 0.5 m/s但不要超过 1 m/s否则会切掉真实的加速早期数据。3.2 拟合主模型参数初值从数据里抠不靠猜拟合使用 scipy 的 curve_fit它内部走的是信赖域反射算法对初值敏感。初值给得越贴近真实收敛越快也越不容易掉进局部最优。以 30 Hz 采样的 100 m 冲刺为例我通常这样设初值参数初值来源典型范围v_max实测速度的 95 分位数612 m/sτ从起跑点到接近 60% vmax 时刻的粗略估计0.54 sγ固定先给 0.0500.5t_p平滑后数据中速度最大值对应的时间29 sσ固定给 3.016 sv_max 用 95 分位数而不是最大值是因为平滑后的数据仍可能有孤立的噪声尖峰最大值容易被单点噪声带偏。τ 的初值可以从「速度达到大约 0.63 v_max 的时间」反推因为单指数上升里 v(τ) ≈ 0.632 v_max。这样设定初值后拟合代码可以写成from scipy.optimize import curve_fit def fit_sprint(t, v): v_max_init np.quantile(v, 0.95) tau_init t[np.argmin(np.abs(v - 0.63 * v_max_init))] t_peak_init t[np.argmax(v)] p0 [v_max_init, tau_init, 0.05, t_peak_init, 3.0] bounds ([5.0, 0.3, 0.0, 0.0, 0.5], [14.0, 5.0, 0.5, 10.0, 8.0]) popt, pcov curve_fit( model_composite, t, v, p0p0, boundsbounds, maxfev40000 ) perr np.sqrt(np.diag(pcov)) return popt, perr popt, perr fit_sprint(t, v) print(v_max {:.3f} ± {:.3f} m/s.format(popt[0], perr[0]))边界约束里 v_max 下限 5 m/s 是针对成年人短跑的下限如果是青少年或普通人群要改成 4 m/s 左右上限 14 m/s 已经高于 100 m 世界级选手的峰值速度留足了余量。γ 的上限给 0.5防止模型把掉速幅度解释得过大而扭曲加速段形状。拟合返回的 pcov 是参数协方差矩阵对角线的平方根就是参数标准误它的数值有没有意义取决于数据是否足以为每个参数提供信息。3.3 拟合失败的三种表现与检查方法拟合输出后第一步不是看参数而是看残差和时间曲线形状。常见的有三种表现提示矩阵接近奇异、参数跑到边界、拟合曲线在起跑段明显穿过数据云下方。矩阵接近奇异通常意味着数据不足以识别全部 5 个参数常见于试跑距离太短只有 60 m或采样率过低。参数跑到边界时比如 τ 顶在 0.3 的下边界多半是起跑段噪声尖峰没有清理干净。起跑段数据不贴合则要考虑给前 0.2 s 的数据降权这个问题会在后面排错章节专门展开。4. 多次试跑的速度特征提取与对比归一化与置信带单次试跑的拟合只能回答「这次跑得怎么样」回答不了「这个运动员的稳定水平是多少」。真实训练场景里同一个运动员会进行多次试跑风速、起跑反应、疲劳程度都会让速度曲线产生差异。这一章处理的是多试跑数据的对比方法。4.1 为什么不能直接对原始速度曲线求平均不同试跑的全程用时不同速度曲线在时间轴上长度不一致直接按时间求平均慢的那次试跑会把曲线向右拉伸快的那次会向左压缩得到的平均曲线既不代表快也不代表慢而是一个两不沾的混合形态。对比多次试跑时先把每条曲线插值到一个公共时间网格上网格范围取所有试跑的最短时长或者取一个固定区间比如 012 st_grid np.linspace(0, min(t_end_list), 500) V np.vstack([ np.interp(t_grid, ti, vi) for ti, vi in zip(t_list, v_list) ]) median_v np.median(V, axis0)公共网格的采样密度取 500 个点在大多数情况下足够再多不会增加信息量因为原始数据是 30 Hz插值到 500 点的时间分辨率已经远高于原始数据。插值方法优先用线性插值PCHIP 在这里反而没必要因为后续还要做非参数置信带估计不需要曲线光滑。4.2 用 bootstrap 给速度曲线加置信带随着试跑次数增多可以用 bootstrap 逐时间点构造置信带。做法是对试跑序号做有放回抽样每次抽样 n 条曲线计算该次抽样的逐点中位数重复 2000 次取每个时间点的 5% 和 95% 分位数作为置信带。这个置信带的含义是「中位数曲线的抽样波动范围」它能直接看出某个时间点的试跑间稳定性rng np.random.default_rng(0) n_boot 2000 n_trials V.shape[0] n_grid V.shape[1] boot_medians np.empty((n_boot, n_grid)) for b in range(n_boot): idx rng.integers(0, n_trials, n_trials) boot_medians[b] np.median(V[idx], axis0) band_low np.quantile(boot_medians, 0.05, axis0) band_high np.quantile(boot_medians, 0.95, axis0)特别注意要先对试跑序号抽样再取中位数而不是对测速点抽样后者会破坏同一试跑内部的速度序列相关性得到的置信带会窄得不真实。2000 次 bootstrap 在数据量小的时候是足够稳定的如果中间章里整个流程跑起来太慢可以降到 1000 次但不要低于这个数否则置信带边界会毛毛糙糙的。4.3 从拟合参数做试跑间对比与异常试跑识别曲线层面的对比之外更常用的是对每次试跑拟合出特征参数然后看参数分布。每次试跑得到一组 v_max、τ、γ、t_p把这些参数汇总后用箱线图或简单的描述统计就能看出稳定性。识别异常试跑的一条实用规则是如果某次试跑的 v_max 低于全部试跑中位数减去 2 倍 MAD或者 t_p 超出实际比赛距离对应的合理窗口就把这次试跑标记为异常并单独复查原始数据。这类异常通常对应起跑失误、途中被干扰或者测速设备丢帧。5. 短跑速度拟合中的 4 个高频坑峰值误判、积分漂移、伪双峰与残差自相关模型选对了拟合代码也跑通了剩下的问题集中在结果解读层面。下面四个坑是我在真实项目里反复踩过的每一条都对应一个具体的数值现象和排查手段。5.1 起跑段反冲毛刺让 τ 估计失真用加权拟合压掉前 0.2 s起跑后 0.10.3 s 内测速数据常常出现一个短暂的速度尖峰来源是脚蹬起跑器后的反弹以及测速仪对突然加速的响应过冲。这个尖峰如果不处理加速时间常数 τ 会被拟合得异常小。处理方式除了滤波外更直接的是在拟合时给前 0.2 s 的数据降权。curve_fit 支持 sigma 参数将权重设为目标函数方差的倒数因此可以把前 0.2 s 的权重视为其他时刻的 5%sigma_w np.ones_like(t) * 0.1 sigma_w[t 0.2] 0.005 # 方差放大 20 倍权重大幅降低 popt, pcov curve_fit(model_composite, t, v, p0p0, sigmasigma_w, boundsbounds, maxfev40000)这里的逻辑不是删除数据而是让拟合算法知道这段数据不可靠。0.2 s 这个阈值对应起跑反冲的持续时间不同运动员差异不大如果测速设备悬挂位置偏低或角度倾斜这个窗口要适当延长到 0.3 s。5.2 积分漂移校验速度积分出来的距离与真实赛程对不上速度曲线拟合完成后一个必须做的验证是把拟合速度曲线对时间积分得到的距离应该和实际赛程高度一致。100 m 比赛拟合出 99 m 或 101 m 很正常但如果偏差超过 2%优先怀疑时间轴的零点不对或者速度单位有问题而不是急着改模型。计算方式如下from scipy.integrate import quad dist_fit, _ quad(lambda tt: model_composite(tt, *popt), 0, t[-1]) dist_error (dist_fit - 100.0) / 100.0 * 100 print(distance error: {:.2f} %.format(dist_error))积分区间从 0 到 t[-1]前提是 t 已经从起跑点归零。如果用的是比赛总用时而不是测速仪的有效数据时长积分区间会偏短误差表现为系统性负偏差。这个检查应该在每次拟合后都执行因为它是发现数据采集问题的低成本手段。5.3 滤波窗口导致的伪双峰峰值时刻被算法拖拽30 Hz 的速度数据如果直接用窗口长度大于 1 s 的 Savitzky-Golay 滤波速度峰值会被明显钝化甚至出现「先达到一个小峰、再跌一点、再升到大峰」的双峰假象。原因很简单窗口太长会把加速段尽头的小波动抹平成两个凸起。处理原则是窗口长度不超过 0.5 s我一般取 0.30.4 s。如果你发现拟合出的 t_p 和从原始数据肉眼看到的峰值时刻相差超过 0.5 s先检查滤波窗口不要怀疑模型。5.4 残差自相关会使参数标准误虚低curve_fit 给出的参数标准误假设残差是独立同分布的但短跑速度数据的残差显然不满足这一点相邻采样点的残差高度相关真实的标准误要比 PCov 给出的数值大。做参数对比时要注意这一点v_max 差 0.1 m/s 未必显著。如果项目需要严格的显著性判断用残差 bootstrap 代替理论标准误做法已经在第 4 章给出只是把重采样对象换成残差序列重新生成拟合数据这样得到的参数分布才是可信的。6. 把短跑速度建模固化成一个可复用脚本流程封装与蒙特卡洛自测前面的章节覆盖了从数据到参数的完整链路但每次分析都重写一遍代码会引入大量隐性错误。最后一步是把整套流程固化成可复用的命令行脚本并且用蒙特卡洛模拟验证拟合器本身在已知真值情况下能还原参数。6.1 封装核心流程为一个 analyze_sprint 函数将预处理、拟合、积分校验和特征提取合并到一个入口函数中输入是原始时间轴和速度序列输出是包含参数、误差和校验结果的字典。以下是精简后的核心实现def analyze_sprint(t, v, distance100.0, fs30.0): # 1. 滤波与起跑点对齐 window int(0.4 * fs) | 1 v savgol_filter(v, window, polyorder3) idx_start np.argmax(v 0.3) t t - t[idx_start] v v[idx_start:] # 2. 加权拟合前 0.2 s 降权 sigma_w np.ones_like(t) * 0.1 sigma_w[t 0.2] 0.005 popt, pcov curve_fit(model_composite, t, v, p0guess_params(t, v), sigmasigma_w, boundsbounds, maxfev40000) # 3. 积分校验 dist_fit, _ quad(lambda tt: model_composite(tt, *popt), 0, t[-1]) # 4. 输出报告 perr np.sqrt(np.diag(pcov)) return { v_max: popt[0], v_max_err: perr[0], tau: popt[1], gamma: popt[2], t_peak: popt[3], dist_error_pct: (dist_fit - distance) / distance * 100 }这里有个细节np.argmax(v 0.3)利用的是布尔数组取 argmax 返回第一个 True 位置这一特性写起来简洁但可读性稍差团队协作时建议改写成显式循环。积分校验虽然已经在前一章出现过但在这个封装函数里它必须是默认步骤而不是可选步骤因为参数输出的物理意义完全取决于拟合曲线是否与实际赛程吻合。6.2 用蒙特卡洛模拟验证拟合器本身可辨识封装的函数是否值得信任需要先在一组已知真值上做测试。做法是生成一条复合模型曲线加上不同水平的白噪声重复拟合数百次统计每个参数估计值的偏倚和标准差。下面是一个最小实现def monte_carlo_check(theta_true, n_sim200, noise0.1): t np.linspace(0, 10, 301) results [] rng np.random.default_rng(42) for _ in range(n_sim): v_true model_composite(t, *theta_true) v_obs v_true rng.normal(0, noise, sizet.shape) try: report analyze_sprint(t, v_obs) results.append(report) except RuntimeError: continue return np.mean([r[v_max] for r in results])运行这段模拟时n_sim 取 200噪声幅度取 0.1 m/s观察 v_max 的均值是否落在真值附近标准差是否合理。如果均值偏离真值超过 0.1 m/s或者超过 20% 的试算拟合失败说明模型在这个数据条件下不可辨识要先回头检查边界约束和初值策略。这个自测流程应该放在项目的最前面作为脚本的一部分而不是分析完数据才开始做因为只有先证明拟合器在合成数据上可靠对真实数据的结论才有依据。本文还有配套的精品资源点击获取 短跑速度曲线拟合：从模型选择到Python实现
简介针对短程赛跑中运动员速度变化问题这份数学建模文档基于Keller赛跑模型建立动态优化模型重点分析生理条件、内外阻力与冲力限制对速度的影响通过微分方程推导速度与路程随时间变化的表达式并配合MATLAB非线性拟合求解最大速度时刻。文档共1个doc文件大小4.69MB包含问题重述、模型假设、公式推导、实际数据对比与结论评价等完整学习资源尤其给出了某届奥运会百米决赛前6名速度数据的理论值与实测值比较误差较小适合数学建模备赛者或运动生物力学初学者参考。该资源已有126人学习文档结构清晰图表与数据表格完整能帮助读者系统理解从微分方程建模、曲线拟合到结果检验的完整流程并可直接借鉴Keller模型的建模思路用于类似运动情景分析。1. 短程赛跑的速度曲线远不是三段论那么整齐短程赛跑的速度曲线不是什么「先加速、后匀速、再减速」的三段论形状。用 30 Hz 以上频响的激光测速仪去测 100 m 冲刺数据点会以一种很不规整的方式抖动起跑后的力反冲、途中跑的步频波动、终点前的小幅掉速全部叠加在一条只有 812 m/s 的曲线上。数学建模要做的事不是把曲线平滑得好看而是用带物理约束的参数方程把「加速快慢」「峰值速度」「掉速幅度」这几个可比较的数字提取出来。这类工作直接服务于运动表现分析比如判断某次试跑是起跑慢了还是后程动力不足。下面按一次完整的数据处理流程展开先选模型再讲拟合再做多次试跑的横向对比最后落到排错与脚本复用。适合接触过 Python 和 numpy、想用曲线拟合处理时序数据的工程人员读。2. 短跑速度-时间模型的 3 个候选方程先决定拟合曲线的「骨架」给运动员速度变化建模第一步不是写拟合代码而是确定曲线的函数形式。曲线方程决定了后面所有参数的物理含义也决定了拟合失败时的排查方向。这里从我做这类数据的经验出发比较三个常用候选单指数上升、指数上升加高斯衰减、非参数样条。2.1 单指数上升起点简单但兜不住后程掉速单指数模型形式上是最简洁的v(t) v_max · (1 - e^(-t / τ))这个方程描述的是「速度从零开始以时间常数 τ 逼近极限速度 v_max」的过程。起跑后的加速段表现得很好但真实短跑在 4070 m 之后速度会从峰值缓慢下降单指数模型缺少这一项拟合出的 v_max 会被后程数据拉低而 τ 会被起跑段主导形成一种两边都不讨好的妥协。所以它更适合做基线模型用来和更复杂的模型做信息准则对比而不是直接作为最终报告用的方程。2.2 指数上升 高斯衰减兼顾加速段与掉速段的实用复合形式实际处理 30 Hz 测速数据时我常用的是一个带衰减项的复合模型v(t) v_max · (1 - e^(-t / τ)) · (1 - γ · e^(-((t - t_p) / σ)^2))这里第一个括号是指数上升部分控制起跑加速的节奏第二个括号是一个反向高斯钟形以 t_p 为中心、σ 为宽度乘上幅度系数 γ 描述后程速度衰减。当 γ0 时退化为单指数模型所以两个模型是嵌套的可以用似然比检验判断衰减项是否显著。这个函数形式并不是哪本标准教材里的命名而是工程上用来同时描述「加速」和「掉速」的复合写法胜在参数非常直观v_max 是理论峰值速度τ 越小加速越快γ 是掉速幅度t_p 是峰值速度出现的时间。2.3 PCHIP 样条不预设形态但参数不可解释如果完全不想预设曲线的数学形态可以用 PCHIP分段三次 Hermite 插值做非参数化拟合。PCHIP 比普通三次样条好在不会在数据跳变处产生过冲速度曲线起跑段那种陡峭上升它也能兜住。但问题在于它的输出是一组插值系数而不是几个可读的参数没办法直接说「这个运动员的加速时间常数是 1.8 s」。用它的正确姿势是作为参考曲线先用 PCHIP 找到峰值速度和峰值时间的近似值再把这些值作为参数模型的初值。2.4 三个候选的取舍与对比模型参数数量可解释性后程掉速拟合稳定性主要用途单指数上升2强不支持高基线对比、起跑段单独分析指数上升 高斯衰减5强支持中全程速度曲线主模型PCHIP 样条随节点数无隐含高峰值位置估计、初值计算这三个模型不是竞争关系而是配合关系。实际流程是先用样条估计峰值位置再用单指数做基线最后用复合模型做正式拟合通过 AIC 确认衰减项的加入是否值得。下面这段代码把三个模型放在同一张图上对比各自的拟合效果方便直观感受差异import numpy as np import matplotlib.pyplot as plt from scipy.interpolate import PCHIP from scipy.optimize import curve_fit def model_single(t, v_max, tau): return v_max * (1 - np.exp(-t / tau)) def model_composite(t, v_max, tau, gamma, t_p, sigma): rise 1 - np.exp(-t / tau) decay 1 - gamma * np.exp(-((t - t_p) / sigma) ** 2) return v_max * rise * decay # 模拟一段 10 秒、30 Hz 的速度数据 fs 30.0 t np.linspace(0, 10, 301) v_true model_composite(t, 11.0, 1.8, 0.08, 6.0, 3.0) rng np.random.default_rng(4) v_obs v_true rng.normal(0, 0.15, sizet.shape) # 单指数拟合 p_single, _ curve_fit(model_single, t, v_obs, p0[10, 2]) # 复合模型拟合 p_comp, _ curve_fit(model_composite, t, v_obs, p0[10, 2, 0.05, 6, 3], maxfev20000) # PCHIP 插值不做平滑作为参考 pchip PCHIP(t, v_obs) t_fine np.linspace(0, 10, 601) plt.figure(figsize(8, 4)) plt.plot(t, v_obs, ., markersize2, alpha0.5, labelobserved) plt.plot(t_fine, model_single(t_fine, *p_single), labelsingle exp) plt.plot(t_fine, model_composite(t_fine, *p_comp), labelexp decay) plt.plot(t_fine, pchip(t_fine), labelPCHIP) plt.xlabel(time (s)); plt.ylabel(velocity (m/s)) plt.legend()代码里复合模型初值中 t_p 给 6 s是因为真实数据的峰值通常出现在 58 s 区间γ 初值给到 0.05 表示先假设轻微掉速。PCHIP 不需要初值它的作用是在参数模型拟合失败时提供一条不坏的可视化曲线帮助你判断是数据问题还是模型问题。3. 用 scipy 拟合运动员速度曲线预处理、初值设定与边界约束模型定下来之后拟合本身反而容易出问题因为测速数据不可能像上一节模拟数据那样干净。真实数据里有起跑的反冲尖峰、有步频带来的周期性波动、还有传感器偶尔丢帧。这一章先处理数据再进入拟合最后看拟合质量如何检验。3.1 速度测量数据的预处理降噪、起点对齐与重采样激光测速仪输出的原始速度序列噪声相当大尤其是起跑后前 2 秒运动员脚掌蹬离起跑器时的反冲会在曲线上留下一个明显的毛刺。直接用原始数据拟合模型会把毛刺当做一个真实的加速特征τ 的估计值会被严重拉偏。我一般先用 Savitzky-Golay 滤波做一次平滑窗口长度取奇数大约是 0.30.5 s 对应的采样点数from scipy.signal import savgol_filter fs 30.0 # 采样率单位 Hz window int(0.4 * fs) | 1 # 保证奇数0.4s 窗口 v_smooth savgol_filter(v_obs, window_lengthwindow, polyorder3) # 判定起跑点速度第一次持续超过 0.3 m/s 的时刻 mask v_smooth 0.3 idx_start np.argmax(mask) t t - t[idx_start] # 时间轴从起跑点重新归零 v v_smooth[idx_start:]这里把时间零点放在起跑瞬间很重要因为复合模型的时间 t 是「起跑后经过的时间」如果时间轴从仪器开机开始算τ 和 t_p 两个参数的物理含义都会错位。0.3 m/s 这个阈值不是硬性标准它的作用是过滤掉传感器零点漂移导致的虚假起始段如果运动员起跑前有准备姿势的微小移动可以把阈值提高到 0.5 m/s但不要超过 1 m/s否则会切掉真实的加速早期数据。3.2 拟合主模型参数初值从数据里抠不靠猜拟合使用 scipy 的 curve_fit它内部走的是信赖域反射算法对初值敏感。初值给得越贴近真实收敛越快也越不容易掉进局部最优。以 30 Hz 采样的 100 m 冲刺为例我通常这样设初值参数初值来源典型范围v_max实测速度的 95 分位数612 m/sτ从起跑点到接近 60% vmax 时刻的粗略估计0.54 sγ固定先给 0.0500.5t_p平滑后数据中速度最大值对应的时间29 sσ固定给 3.016 sv_max 用 95 分位数而不是最大值是因为平滑后的数据仍可能有孤立的噪声尖峰最大值容易被单点噪声带偏。τ 的初值可以从「速度达到大约 0.63 v_max 的时间」反推因为单指数上升里 v(τ) ≈ 0.632 v_max。这样设定初值后拟合代码可以写成from scipy.optimize import curve_fit def fit_sprint(t, v): v_max_init np.quantile(v, 0.95) tau_init t[np.argmin(np.abs(v - 0.63 * v_max_init))] t_peak_init t[np.argmax(v)] p0 [v_max_init, tau_init, 0.05, t_peak_init, 3.0] bounds ([5.0, 0.3, 0.0, 0.0, 0.5], [14.0, 5.0, 0.5, 10.0, 8.0]) popt, pcov curve_fit( model_composite, t, v, p0p0, boundsbounds, maxfev40000 ) perr np.sqrt(np.diag(pcov)) return popt, perr popt, perr fit_sprint(t, v) print(v_max {:.3f} ± {:.3f} m/s.format(popt[0], perr[0]))边界约束里 v_max 下限 5 m/s 是针对成年人短跑的下限如果是青少年或普通人群要改成 4 m/s 左右上限 14 m/s 已经高于 100 m 世界级选手的峰值速度留足了余量。γ 的上限给 0.5防止模型把掉速幅度解释得过大而扭曲加速段形状。拟合返回的 pcov 是参数协方差矩阵对角线的平方根就是参数标准误它的数值有没有意义取决于数据是否足以为每个参数提供信息。3.3 拟合失败的三种表现与检查方法拟合输出后第一步不是看参数而是看残差和时间曲线形状。常见的有三种表现提示矩阵接近奇异、参数跑到边界、拟合曲线在起跑段明显穿过数据云下方。矩阵接近奇异通常意味着数据不足以识别全部 5 个参数常见于试跑距离太短只有 60 m或采样率过低。参数跑到边界时比如 τ 顶在 0.3 的下边界多半是起跑段噪声尖峰没有清理干净。起跑段数据不贴合则要考虑给前 0.2 s 的数据降权这个问题会在后面排错章节专门展开。4. 多次试跑的速度特征提取与对比归一化与置信带单次试跑的拟合只能回答「这次跑得怎么样」回答不了「这个运动员的稳定水平是多少」。真实训练场景里同一个运动员会进行多次试跑风速、起跑反应、疲劳程度都会让速度曲线产生差异。这一章处理的是多试跑数据的对比方法。4.1 为什么不能直接对原始速度曲线求平均不同试跑的全程用时不同速度曲线在时间轴上长度不一致直接按时间求平均慢的那次试跑会把曲线向右拉伸快的那次会向左压缩得到的平均曲线既不代表快也不代表慢而是一个两不沾的混合形态。对比多次试跑时先把每条曲线插值到一个公共时间网格上网格范围取所有试跑的最短时长或者取一个固定区间比如 012 st_grid np.linspace(0, min(t_end_list), 500) V np.vstack([ np.interp(t_grid, ti, vi) for ti, vi in zip(t_list, v_list) ]) median_v np.median(V, axis0)公共网格的采样密度取 500 个点在大多数情况下足够再多不会增加信息量因为原始数据是 30 Hz插值到 500 点的时间分辨率已经远高于原始数据。插值方法优先用线性插值PCHIP 在这里反而没必要因为后续还要做非参数置信带估计不需要曲线光滑。4.2 用 bootstrap 给速度曲线加置信带随着试跑次数增多可以用 bootstrap 逐时间点构造置信带。做法是对试跑序号做有放回抽样每次抽样 n 条曲线计算该次抽样的逐点中位数重复 2000 次取每个时间点的 5% 和 95% 分位数作为置信带。这个置信带的含义是「中位数曲线的抽样波动范围」它能直接看出某个时间点的试跑间稳定性rng np.random.default_rng(0) n_boot 2000 n_trials V.shape[0] n_grid V.shape[1] boot_medians np.empty((n_boot, n_grid)) for b in range(n_boot): idx rng.integers(0, n_trials, n_trials) boot_medians[b] np.median(V[idx], axis0) band_low np.quantile(boot_medians, 0.05, axis0) band_high np.quantile(boot_medians, 0.95, axis0)特别注意要先对试跑序号抽样再取中位数而不是对测速点抽样后者会破坏同一试跑内部的速度序列相关性得到的置信带会窄得不真实。2000 次 bootstrap 在数据量小的时候是足够稳定的如果中间章里整个流程跑起来太慢可以降到 1000 次但不要低于这个数否则置信带边界会毛毛糙糙的。4.3 从拟合参数做试跑间对比与异常试跑识别曲线层面的对比之外更常用的是对每次试跑拟合出特征参数然后看参数分布。每次试跑得到一组 v_max、τ、γ、t_p把这些参数汇总后用箱线图或简单的描述统计就能看出稳定性。识别异常试跑的一条实用规则是如果某次试跑的 v_max 低于全部试跑中位数减去 2 倍 MAD或者 t_p 超出实际比赛距离对应的合理窗口就把这次试跑标记为异常并单独复查原始数据。这类异常通常对应起跑失误、途中被干扰或者测速设备丢帧。5. 短跑速度拟合中的 4 个高频坑峰值误判、积分漂移、伪双峰与残差自相关模型选对了拟合代码也跑通了剩下的问题集中在结果解读层面。下面四个坑是我在真实项目里反复踩过的每一条都对应一个具体的数值现象和排查手段。5.1 起跑段反冲毛刺让 τ 估计失真用加权拟合压掉前 0.2 s起跑后 0.10.3 s 内测速数据常常出现一个短暂的速度尖峰来源是脚蹬起跑器后的反弹以及测速仪对突然加速的响应过冲。这个尖峰如果不处理加速时间常数 τ 会被拟合得异常小。处理方式除了滤波外更直接的是在拟合时给前 0.2 s 的数据降权。curve_fit 支持 sigma 参数将权重设为目标函数方差的倒数因此可以把前 0.2 s 的权重视为其他时刻的 5%sigma_w np.ones_like(t) * 0.1 sigma_w[t 0.2] 0.005 # 方差放大 20 倍权重大幅降低 popt, pcov curve_fit(model_composite, t, v, p0p0, sigmasigma_w, boundsbounds, maxfev40000)这里的逻辑不是删除数据而是让拟合算法知道这段数据不可靠。0.2 s 这个阈值对应起跑反冲的持续时间不同运动员差异不大如果测速设备悬挂位置偏低或角度倾斜这个窗口要适当延长到 0.3 s。5.2 积分漂移校验速度积分出来的距离与真实赛程对不上速度曲线拟合完成后一个必须做的验证是把拟合速度曲线对时间积分得到的距离应该和实际赛程高度一致。100 m 比赛拟合出 99 m 或 101 m 很正常但如果偏差超过 2%优先怀疑时间轴的零点不对或者速度单位有问题而不是急着改模型。计算方式如下from scipy.integrate import quad dist_fit, _ quad(lambda tt: model_composite(tt, *popt), 0, t[-1]) dist_error (dist_fit - 100.0) / 100.0 * 100 print(distance error: {:.2f} %.format(dist_error))积分区间从 0 到 t[-1]前提是 t 已经从起跑点归零。如果用的是比赛总用时而不是测速仪的有效数据时长积分区间会偏短误差表现为系统性负偏差。这个检查应该在每次拟合后都执行因为它是发现数据采集问题的低成本手段。5.3 滤波窗口导致的伪双峰峰值时刻被算法拖拽30 Hz 的速度数据如果直接用窗口长度大于 1 s 的 Savitzky-Golay 滤波速度峰值会被明显钝化甚至出现「先达到一个小峰、再跌一点、再升到大峰」的双峰假象。原因很简单窗口太长会把加速段尽头的小波动抹平成两个凸起。处理原则是窗口长度不超过 0.5 s我一般取 0.30.4 s。如果你发现拟合出的 t_p 和从原始数据肉眼看到的峰值时刻相差超过 0.5 s先检查滤波窗口不要怀疑模型。5.4 残差自相关会使参数标准误虚低curve_fit 给出的参数标准误假设残差是独立同分布的但短跑速度数据的残差显然不满足这一点相邻采样点的残差高度相关真实的标准误要比 PCov 给出的数值大。做参数对比时要注意这一点v_max 差 0.1 m/s 未必显著。如果项目需要严格的显著性判断用残差 bootstrap 代替理论标准误做法已经在第 4 章给出只是把重采样对象换成残差序列重新生成拟合数据这样得到的参数分布才是可信的。6. 把短跑速度建模固化成一个可复用脚本流程封装与蒙特卡洛自测前面的章节覆盖了从数据到参数的完整链路但每次分析都重写一遍代码会引入大量隐性错误。最后一步是把整套流程固化成可复用的命令行脚本并且用蒙特卡洛模拟验证拟合器本身在已知真值情况下能还原参数。6.1 封装核心流程为一个 analyze_sprint 函数将预处理、拟合、积分校验和特征提取合并到一个入口函数中输入是原始时间轴和速度序列输出是包含参数、误差和校验结果的字典。以下是精简后的核心实现def analyze_sprint(t, v, distance100.0, fs30.0): # 1. 滤波与起跑点对齐 window int(0.4 * fs) | 1 v savgol_filter(v, window, polyorder3) idx_start np.argmax(v 0.3) t t - t[idx_start] v v[idx_start:] # 2. 加权拟合前 0.2 s 降权 sigma_w np.ones_like(t) * 0.1 sigma_w[t 0.2] 0.005 popt, pcov curve_fit(model_composite, t, v, p0guess_params(t, v), sigmasigma_w, boundsbounds, maxfev40000) # 3. 积分校验 dist_fit, _ quad(lambda tt: model_composite(tt, *popt), 0, t[-1]) # 4. 输出报告 perr np.sqrt(np.diag(pcov)) return { v_max: popt[0], v_max_err: perr[0], tau: popt[1], gamma: popt[2], t_peak: popt[3], dist_error_pct: (dist_fit - distance) / distance * 100 }这里有个细节np.argmax(v 0.3)利用的是布尔数组取 argmax 返回第一个 True 位置这一特性写起来简洁但可读性稍差团队协作时建议改写成显式循环。积分校验虽然已经在前一章出现过但在这个封装函数里它必须是默认步骤而不是可选步骤因为参数输出的物理意义完全取决于拟合曲线是否与实际赛程吻合。6.2 用蒙特卡洛模拟验证拟合器本身可辨识封装的函数是否值得信任需要先在一组已知真值上做测试。做法是生成一条复合模型曲线加上不同水平的白噪声重复拟合数百次统计每个参数估计值的偏倚和标准差。下面是一个最小实现def monte_carlo_check(theta_true, n_sim200, noise0.1): t np.linspace(0, 10, 301) results [] rng np.random.default_rng(42) for _ in range(n_sim): v_true model_composite(t, *theta_true) v_obs v_true rng.normal(0, noise, sizet.shape) try: report analyze_sprint(t, v_obs) results.append(report) except RuntimeError: continue return np.mean([r[v_max] for r in results])运行这段模拟时n_sim 取 200噪声幅度取 0.1 m/s观察 v_max 的均值是否落在真值附近标准差是否合理。如果均值偏离真值超过 0.1 m/s或者超过 20% 的试算拟合失败说明模型在这个数据条件下不可辨识要先回头检查边界约束和初值策略。这个自测流程应该放在项目的最前面作为脚本的一部分而不是分析完数据才开始做因为只有先证明拟合器在合成数据上可靠对真实数据的结论才有依据。本文还有配套的精品资源点击获取