尧图网页设计 
简介本资源是2024年美国大学生数学建模竞赛MCM/ICMProblem E的完整参赛解决方案聚焦财产保险可持续性评估这一现实金融议题面向数学建模初学者、保险科技学习者及Python数据分析进阶用户适用于课程设计、毕业设计与竞赛备赛。压缩包共43个文件含7个核心Python脚本涵盖灰色预测、模糊综合评价、支持向量机建模等、16张结果可视化图表PNG/SVG格式、3个关键数据表XLSX、2份PDF文档含赛题原文与技术报告、1个模型参数文件PKL及LaTeX源码TEX整体18.43MB结构清晰、模块可拆解复用。已有110人下载学习提供从数据预处理、多模型构建GM(1,1)、SVM、PCA-AUC分析、灵敏度检验到可视化呈现的全流程实现附带README说明与分步执行逻辑便于理解建模思路、复现关键结果并拓展至其他风险评估场景。1. 财产保险可持续性不是财务报表里的空话而是能用 Python 算出来的动态平衡2024年美国大学生数学建模竞赛MCMC题聚焦财产保险的可持续性——这不是在讨论公司要不要发CSR报告而是直击核心当一场区域性暴雨导致某地37%的住宅屋顶受损、车险报案量单日激增210%、再保险合约触发自动摊赔条款时这家保险公司未来18个月的偿付能力充足率是否仍高于150%能否在不提价20%、不拒保高风险区域的前提下维持承保利润为正这类问题无法靠Excel滚动预测解决必须构建包含损失分布建模、资本缓冲动态计算、费率弹性响应、再保险结构嵌套的多层系统。Python 成为此类建模的首选并非因为语法简单而是其 SciPy 的极值分布拟合、Statsmodels 的广义线性混合模型GLMM、PyMC 的贝叶斯分层建模、以及 NetworkX 对再保险合约网络的拓扑分析能力恰好覆盖了财产保险可持续性的四个刚性技术断点。本文面向已掌握 Pandas 基础、正参与数模集训或准备 MCM/ICM 的本科生不讲“什么是保险”只拆解如何用 67 行核心代码跑通从历史赔案数据到监管资本阈值预警的完整链路。2. 用 Python 构建财产保险可持续性评估的最小可行模型财产保险可持续性建模的本质是建立“损失生成—资本消耗—费率调节—再保缓冲”四环节的闭环反馈系统。常见误区是直接套用 Logistic 回归预测“是否破产”这忽略了保险经营中时间维度的路径依赖性与监管规则的硬约束。我们采用分层建模策略底层用极值理论拟合巨灾损失尾部中层用随机过程模拟资本充足率演化上层嵌入监管阈值触发的费率调整规则。这种结构既满足 MCM 对“可解释性”的硬性要求又具备向更复杂模型如加入气候模型输出作为协变量扩展的接口。2.1 从原始赔案数据中提取可建模的损失特征MCM 官方数据包通常提供结构化赔案表claim_id, policy_id, loss_amount, claim_date, zip_code, construction_type但原始字段无法直接输入模型。关键预处理步骤有三第一剔除明显异常值——不是简单用 3σ 法则而是采用Hampel 标识器它对时间序列中的脉冲噪声更鲁棒import numpy as np from scipy import signal def hampel_filter(data, window_size5, n_sigmas3): 对损失金额序列进行Hampel滤波window_size为滑动窗口长度 median signal.medfilt(data, kernel_sizewindow_size) mad np.median(np.abs(data - median)) threshold n_sigmas * 1.4826 * mad # MAD转标准差的常数 outliers np.abs(data - median) threshold data_clean data.copy() data_clean[outliers] median[outliers] return data_clean # 应用于某地区2020-2023年季度累计损失额序列 quarterly_losses np.array([12.5, 14.2, 13.8, 98.7, 15.1, 16.3, 14.9, 152.4]) # 单位百万美元 cleaned hampel_filter(quarterly_losses) print(f原始序列: {quarterly_losses}) print(f清洗后: {cleaned}) # 输出中98.7和152.4被修正为邻近中位数提示signal.medfilt要求window_size为奇数且必须大于等于31.4826是将中位绝对偏差MAD转换为正态分布标准差的理论系数不可省略。第二构造损失频率与严重度的联合特征。财产险需分离“发生多少次事故”频率和“每次赔多少钱”严重度。我们按保单年度聚合import pandas as pd # 假设df_claims为原始赔案DataFrame df_claims[claim_year] pd.to_datetime(df_claims[claim_date]).dt.year df_policy_annual df_claims.groupby([policy_id, claim_year]).agg( claim_count(claim_id, count), total_loss(loss_amount, sum) ).reset_index() # 计算每份保单每年的平均损失严重度避免零索赔保单干扰 df_policy_annual[severity] np.where( df_policy_annual[claim_count] 0, df_policy_annual[total_loss] / df_policy_annual[claim_count], np.nan )第三引入地理风险因子。ZIP Code 本身是分类变量需转换为可建模的连续指标。我们采用风险暴露加权法# 假设risk_map为外部加载的地理风险指数表zip_code - flood_risk_score, wind_risk_score df_merged df_policy_annual.merge(risk_map, onzip_code, howleft) # 构造复合风险指数flood权重0.6wind权重0.4经Z-score标准化 df_merged[risk_index] ( 0.6 * (df_merged[flood_risk_score] - df_merged[flood_risk_score].mean()) / df_merged[flood_risk_score].std() 0.4 * (df_merged[wind_risk_score] - df_merged[wind_risk_score].mean()) / df_merged[wind_risk_score].std() )2.2 用广义帕累托分布GPD拟合损失尾部并计算VaR财产险巨灾损失具有尖峰厚尾特性正态分布完全失效。GPD 是极值理论中描述超额损失的标准工具其累积分布函数为$$F(x) 1 - \left(1 \xi \frac{x-u}{\sigma}\right)^{-1/\xi}, \quad x u$$其中 $u$ 是阈值$\sigma 0$ 是尺度参数$\xi$ 是形状参数$\xi 0$ 表示重尾$\xi 0$ 退化为指数分布。MCM 评审关注你是否合理选择 $u$而非仅调用.fit()。from scipy.stats import genpareto import matplotlib.pyplot as plt # 步骤1确定阈值u —— 使用平均超额图Mean Excess Plot losses df_policy_annual[total_loss].dropna().values thresholds np.linspace(np.percentile(losses, 90), np.percentile(losses, 99), 20) excess_means [np.mean(losses[losses t] - t) for t in thresholds] plt.figure(figsize(8, 4)) plt.plot(thresholds, excess_means, o-) plt.xlabel(Threshold u (million USD)) plt.ylabel(Mean Excess over u) plt.title(Mean Excess Plot for Loss Tail) plt.grid(True) plt.show() # 观察图中直线段起始点取u 25.0即92.5%分位数 u_threshold np.percentile(losses, 92.5) # 步骤2拟合GPD并计算99.5%置信水平下的VaR监管常用阈值 excesses losses[losses u_threshold] - u_threshold shape, loc, scale genpareto.fit(excesses, floc0) # 强制loc0符合GPD定义 # VaR_99.5% u (scale/shape) * [(1-0.995)^(-shape) - 1] var_995 u_threshold (scale / shape) * ((1 - 0.995)**(-shape) - 1) print(fGPD拟合结果: shape{shape:.4f}, scale{scale:.4f}) print(f99.5% VaR (百万美元): {var_995:.2f})注意genpareto.fit()的floc0参数至关重要它确保位置参数固定为0否则拟合可能发散shape若为负值说明尾部有界此时应改用其他分布如Beta。2.3 构建资本充足率动态演化模型可持续性的核心指标是资本充足率CAR 可用资本 / 最低资本要求。我们用随机差分方程模拟其季度演化$$CAR_{t1} CAR_t \alpha \cdot (Premium_t - Loss_t - Expenses_t) / Capital_t - \beta \cdot \Delta Reinsurance_t$$其中 $\alpha$ 是资本转化效率系数通常0.8~0.95$\beta$ 是再保险成本敏感度0.1~0.3。该模型比静态比率更有说服力因为它揭示了“一次大灾后需要几个季度恢复”。def simulate_car_evolution(initial_car2.1, n_quarters24, premium_growth0.015, loss_volatility0.25): 模拟资本充足率24个季度演化 initial_car: 初始CAR值如2.1表示210% premium_growth: 季度保费增长率小数如0.0151.5% loss_volatility: 损失波动率控制GPD抽样离散度 car_history [initial_car] capital_base 1000 # 初始资本基数百万美元 for t in range(n_quarters): # 保费收入按增长率递增叠加±5%随机扰动 premium 100 * (1 premium_growth)**t * (1 np.random.uniform(-0.05, 0.05)) # 损失支出从GPD抽样但限制单季损失不超过资本的30% loss_sample genpareto.rvs(shape, loc0, scalescale, size1)[0] u_threshold loss min(loss_sample, 0.3 * capital_base) # 防止单季清零 # 运营费用保费的22% expenses 0.22 * premium # 再保险变动当CAR1.8时自动增加再保分出成本上升15% reinsurance_delta 0.0 if car_history[-1] 1.8 else 0.15 # 更新CAR公式中α0.92, β0.25 delta_car 0.92 * (premium - loss - expenses) / capital_base - 0.25 * reinsurance_delta new_car car_history[-1] delta_car car_history.append(max(new_car, 0.5)) # CAR不低于50%否则视为技术性违约 return np.array(car_history) # 运行100次蒙特卡洛模拟 car_simulations np.array([simulate_car_evolution() for _ in range(100)]) car_mean np.mean(car_simulations, axis0) car_5pct np.percentile(car_simulations, 5, axis0) # 5%分位数监管关注底线 plt.figure(figsize(10, 5)) plt.plot(car_mean, b-, label平均CAR) plt.fill_between(range(25), car_5pct, alpha0.3, colorred, label5%分位数) plt.axhline(y1.5, colork, linestyle--, label监管阈值150%) plt.xlabel(Quarter) plt.ylabel(Capital Adequacy Ratio) plt.legend() plt.title(Capital Adequacy Ratio Evolution (100 Simulations)) plt.grid(True) plt.show()3. 将可持续性转化为可执行的费率优化策略建模的终点不是画一条CAR曲线而是回答“下季度对沿海地区砖混结构住宅基准费率应上调多少才能使CAR_5pct稳定在1.5以上” 这要求将统计模型与精算定价逻辑打通。我们采用梯度提升树XGBoost作为连接层输入是风险因子地理指数、建筑类型、免赔额、市场因子同业费率、通胀率输出是费率调整系数。关键创新在于损失函数不是常规的MSE而是CAR约束下的定制损失。3.1 构造带监管约束的XGBoost训练目标传统XGBoost回归会最小化预测误差但保险定价需满足调整后费率必须使模拟CAR ≥ 1.5。我们设计两阶段目标第一阶段用历史数据训练基础费率模型第二阶段在基础模型上叠加一个“约束校正项”该修正项由另一个XGBoost学习其标签是若当前费率导致CAR 1.5则标签为正需提价否则为0。import xgboost as xgb from sklearn.model_selection import train_test_split # 假设df_pricing为定价建模数据集含features列和target列历史费率 X df_pricing[feature_columns] y df_pricing[rate_per_thousand] # 第一阶段基础费率模型 X_train, X_test, y_train, y_test train_test_split(X, y, test_size0.2, random_state42) base_model xgb.XGBRegressor(n_estimators100, max_depth4, learning_rate0.1) base_model.fit(X_train, y_train) # 第二阶段约束校正模型 —— 标签为二元1需提价以满足CAR≥1.50无需 # 关键标签生成需调用2.3节的CAR模拟函数 def generate_constraint_label(row): # 用row中的风险因子模拟CAR若5%分位CAR1.5则返回1 simulated_car simulate_car_evolution( initial_car2.0, premium_growth0.015, loss_volatility0.25 ) return 1 if np.percentile(simulated_car, 5) 1.5 else 0 # 为训练集每行生成标签实际中需批量计算此处简化 y_constraint [generate_constraint_label(row) for _, row in X_train.iterrows()] constraint_model xgb.XGBClassifier(n_estimators50, max_depth3, learning_rate0.15) constraint_model.fit(X_train, y_constraint) # 预测时基础费率 × (1 0.15 × constraint_prediction) y_pred_base base_model.predict(X_test) y_pred_constraint constraint_model.predict(X_test) final_rate y_pred_base * (1 0.15 * y_pred_constraint)3.2 可视化费率调整的地理热力图与可持续性收益MCM 论文要求空间可视化。我们用 GeoPandas 绘制 ZIP Code 级别费率调整幅度并叠加 CAR 改善效果import geopandas as gpd import contextily as ctx # 加载美国ZIP Code地理边界需提前下载shp文件 gdf_zip gpd.read_file(tl_2023_us_zcta520.shp) gdf_zip gdf_zip.to_crs(epsg3857) # 转为Web Mercator投影 # 合并费率调整数据假设df_rates含zip_code和rate_change_pct列 gdf_merged gdf_zip.merge(df_rates, left_onZCTA5CE20, right_onzip_code, howleft) fig, ax plt.subplots(1, 2, figsize(18, 8)) # 左图费率调整热力图 gdf_merged.plot(columnrate_change_pct, cmapRdYlBu_r, legendTrue, axax[0], missing_kwds{color: lightgrey}) ax[0].set_title(Rate Adjustment by ZIP Code (%)) ctx.add_basemap(ax[0], sourcectx.providers.OpenStreetMap.Mapnik) # 右图CAR改善对比调整前vs调整后 car_before car_simulations[0] # 任选一次模拟 car_after simulate_car_evolution( initial_car2.0, premium_growth0.015 0.003, # 假设整体费率提升0.3% loss_volatility0.25 ) ax[1].plot(car_before, r--, labelCAR before adjustment) ax[1].plot(car_after, g-, labelCAR after adjustment) ax[1].axhline(y1.5, colork, linestyle:, labelRegulatory threshold) ax[1].set_xlabel(Quarter) ax[1].set_ylabel(Capital Adequacy Ratio) ax[1].legend() ax[1].grid(True) ax[1].set_title(CAR Evolution: Impact of Rate Adjustment) plt.tight_layout() plt.show()提示contextily的 basemap 需联网加载若离线环境可替换为gdf_zip.boundary.plot(axax[0])仅绘制行政边界。4. 在MCM限时场景下快速验证模型可靠性的3个硬核技巧竞赛中没有时间做全量交叉验证。我们用三个可在30分钟内完成的实证检验直接回应评审最可能质疑的点模型是否过拟合参数是否合理结论是否稳健4.1 用“反事实删除法”检验关键参数敏感性评审常问“如果GPD的shape参数变化±0.1VaR会偏移多少” 我们不重新拟合100次而是用解析法快速计算# GPD的VaR对shape的导数解析解 def var_sensitivity_to_shape(shape, scale, u, p0.995): 计算VaR_p对shape的偏导数 term (1 - p)**(-shape) d_term_d_shape -term * np.log(1 - p) return (scale / shape**2) * (term - 1) - (scale / shape) * d_term_d_shape sens var_sensitivity_to_shape(shape, scale, u_threshold) print(fVaR_99.5%对shape的敏感度: {sens:.2f} (每单位shape变化导致VaR变化)) # 若sens 50说明VaR对shape极度敏感需在论文中讨论shape的贝叶斯估计4.2 用Bootstrap重采样验证CAR阈值的统计显著性判断“CAR_5pct是否真的≥1.5”不能只看点估计。我们用Bootstrap计算其95%置信区间from sklearn.utils import resample def bootstrap_car_ci(simulations, alpha0.05): 对CAR序列的5%分位数计算置信区间 n_sim simulations.shape[0] n_boot 1000 boot_5pct [] for _ in range(n_boot): # 从100次模拟中重采样100次允许重复 indices np.random.randint(0, n_sim, n_sim) boot_sim simulations[indices] boot_5pct.append(np.percentile(boot_sim[:, -1], 5)) # 取最终季度的5%分位 lower np.percentile(boot_5pct, alpha/2 * 100) upper np.percentile(boot_5pct, (1 - alpha/2) * 100) return lower, upper lower_ci, upper_ci bootstrap_car_ci(car_simulations) print(fCAR_5pct at quarter 24: [{lower_ci:.3f}, {upper_ci:.3f}]) # 若区间完全在1.5上方可写“在95%置信水平下可持续性成立”4.3 构建“压力测试仪表盘”一键生成关键图表把上述所有检验封装成函数输入新数据即可输出PDF报告def generate_stress_dashboard(new_data_path): 输入新赔案CSV输出含4张图的PDF df_new pd.read_csv(new_data_path) # ... 执行2.1~2.3节全部流程 ... # 生成四图1. Mean Excess Plot 2. CAR模拟分布 3. 费率热力图 4. 敏感性矩阵 fig, axes plt.subplots(2, 2, figsize(16, 12)) # ... 绘图代码 ... plt.savefig(stress_test_report.pdf, bbox_inchestight) print(压力测试报告已生成: stress_test_report.pdf) # 竞赛中只需一行命令 # generate_stress_dashboard(data/2024Q1_claims.csv)注意bbox_inchestight可自动裁剪图表白边避免PDF中出现大片空白这是MCM排版的隐形加分项。本文还有配套的精品资源点击获取 Python构建财产保险可持续性动态评估模型
简介本资源是2024年美国大学生数学建模竞赛MCM/ICMProblem E的完整参赛解决方案聚焦财产保险可持续性评估这一现实金融议题面向数学建模初学者、保险科技学习者及Python数据分析进阶用户适用于课程设计、毕业设计与竞赛备赛。压缩包共43个文件含7个核心Python脚本涵盖灰色预测、模糊综合评价、支持向量机建模等、16张结果可视化图表PNG/SVG格式、3个关键数据表XLSX、2份PDF文档含赛题原文与技术报告、1个模型参数文件PKL及LaTeX源码TEX整体18.43MB结构清晰、模块可拆解复用。已有110人下载学习提供从数据预处理、多模型构建GM(1,1)、SVM、PCA-AUC分析、灵敏度检验到可视化呈现的全流程实现附带README说明与分步执行逻辑便于理解建模思路、复现关键结果并拓展至其他风险评估场景。1. 财产保险可持续性不是财务报表里的空话而是能用 Python 算出来的动态平衡2024年美国大学生数学建模竞赛MCMC题聚焦财产保险的可持续性——这不是在讨论公司要不要发CSR报告而是直击核心当一场区域性暴雨导致某地37%的住宅屋顶受损、车险报案量单日激增210%、再保险合约触发自动摊赔条款时这家保险公司未来18个月的偿付能力充足率是否仍高于150%能否在不提价20%、不拒保高风险区域的前提下维持承保利润为正这类问题无法靠Excel滚动预测解决必须构建包含损失分布建模、资本缓冲动态计算、费率弹性响应、再保险结构嵌套的多层系统。Python 成为此类建模的首选并非因为语法简单而是其 SciPy 的极值分布拟合、Statsmodels 的广义线性混合模型GLMM、PyMC 的贝叶斯分层建模、以及 NetworkX 对再保险合约网络的拓扑分析能力恰好覆盖了财产保险可持续性的四个刚性技术断点。本文面向已掌握 Pandas 基础、正参与数模集训或准备 MCM/ICM 的本科生不讲“什么是保险”只拆解如何用 67 行核心代码跑通从历史赔案数据到监管资本阈值预警的完整链路。2. 用 Python 构建财产保险可持续性评估的最小可行模型财产保险可持续性建模的本质是建立“损失生成—资本消耗—费率调节—再保缓冲”四环节的闭环反馈系统。常见误区是直接套用 Logistic 回归预测“是否破产”这忽略了保险经营中时间维度的路径依赖性与监管规则的硬约束。我们采用分层建模策略底层用极值理论拟合巨灾损失尾部中层用随机过程模拟资本充足率演化上层嵌入监管阈值触发的费率调整规则。这种结构既满足 MCM 对“可解释性”的硬性要求又具备向更复杂模型如加入气候模型输出作为协变量扩展的接口。2.1 从原始赔案数据中提取可建模的损失特征MCM 官方数据包通常提供结构化赔案表claim_id, policy_id, loss_amount, claim_date, zip_code, construction_type但原始字段无法直接输入模型。关键预处理步骤有三第一剔除明显异常值——不是简单用 3σ 法则而是采用Hampel 标识器它对时间序列中的脉冲噪声更鲁棒import numpy as np from scipy import signal def hampel_filter(data, window_size5, n_sigmas3): 对损失金额序列进行Hampel滤波window_size为滑动窗口长度 median signal.medfilt(data, kernel_sizewindow_size) mad np.median(np.abs(data - median)) threshold n_sigmas * 1.4826 * mad # MAD转标准差的常数 outliers np.abs(data - median) threshold data_clean data.copy() data_clean[outliers] median[outliers] return data_clean # 应用于某地区2020-2023年季度累计损失额序列 quarterly_losses np.array([12.5, 14.2, 13.8, 98.7, 15.1, 16.3, 14.9, 152.4]) # 单位百万美元 cleaned hampel_filter(quarterly_losses) print(f原始序列: {quarterly_losses}) print(f清洗后: {cleaned}) # 输出中98.7和152.4被修正为邻近中位数提示signal.medfilt要求window_size为奇数且必须大于等于31.4826是将中位绝对偏差MAD转换为正态分布标准差的理论系数不可省略。第二构造损失频率与严重度的联合特征。财产险需分离“发生多少次事故”频率和“每次赔多少钱”严重度。我们按保单年度聚合import pandas as pd # 假设df_claims为原始赔案DataFrame df_claims[claim_year] pd.to_datetime(df_claims[claim_date]).dt.year df_policy_annual df_claims.groupby([policy_id, claim_year]).agg( claim_count(claim_id, count), total_loss(loss_amount, sum) ).reset_index() # 计算每份保单每年的平均损失严重度避免零索赔保单干扰 df_policy_annual[severity] np.where( df_policy_annual[claim_count] 0, df_policy_annual[total_loss] / df_policy_annual[claim_count], np.nan )第三引入地理风险因子。ZIP Code 本身是分类变量需转换为可建模的连续指标。我们采用风险暴露加权法# 假设risk_map为外部加载的地理风险指数表zip_code - flood_risk_score, wind_risk_score df_merged df_policy_annual.merge(risk_map, onzip_code, howleft) # 构造复合风险指数flood权重0.6wind权重0.4经Z-score标准化 df_merged[risk_index] ( 0.6 * (df_merged[flood_risk_score] - df_merged[flood_risk_score].mean()) / df_merged[flood_risk_score].std() 0.4 * (df_merged[wind_risk_score] - df_merged[wind_risk_score].mean()) / df_merged[wind_risk_score].std() )2.2 用广义帕累托分布GPD拟合损失尾部并计算VaR财产险巨灾损失具有尖峰厚尾特性正态分布完全失效。GPD 是极值理论中描述超额损失的标准工具其累积分布函数为$$F(x) 1 - \left(1 \xi \frac{x-u}{\sigma}\right)^{-1/\xi}, \quad x u$$其中 $u$ 是阈值$\sigma 0$ 是尺度参数$\xi$ 是形状参数$\xi 0$ 表示重尾$\xi 0$ 退化为指数分布。MCM 评审关注你是否合理选择 $u$而非仅调用.fit()。from scipy.stats import genpareto import matplotlib.pyplot as plt # 步骤1确定阈值u —— 使用平均超额图Mean Excess Plot losses df_policy_annual[total_loss].dropna().values thresholds np.linspace(np.percentile(losses, 90), np.percentile(losses, 99), 20) excess_means [np.mean(losses[losses t] - t) for t in thresholds] plt.figure(figsize(8, 4)) plt.plot(thresholds, excess_means, o-) plt.xlabel(Threshold u (million USD)) plt.ylabel(Mean Excess over u) plt.title(Mean Excess Plot for Loss Tail) plt.grid(True) plt.show() # 观察图中直线段起始点取u 25.0即92.5%分位数 u_threshold np.percentile(losses, 92.5) # 步骤2拟合GPD并计算99.5%置信水平下的VaR监管常用阈值 excesses losses[losses u_threshold] - u_threshold shape, loc, scale genpareto.fit(excesses, floc0) # 强制loc0符合GPD定义 # VaR_99.5% u (scale/shape) * [(1-0.995)^(-shape) - 1] var_995 u_threshold (scale / shape) * ((1 - 0.995)**(-shape) - 1) print(fGPD拟合结果: shape{shape:.4f}, scale{scale:.4f}) print(f99.5% VaR (百万美元): {var_995:.2f})注意genpareto.fit()的floc0参数至关重要它确保位置参数固定为0否则拟合可能发散shape若为负值说明尾部有界此时应改用其他分布如Beta。2.3 构建资本充足率动态演化模型可持续性的核心指标是资本充足率CAR 可用资本 / 最低资本要求。我们用随机差分方程模拟其季度演化$$CAR_{t1} CAR_t \alpha \cdot (Premium_t - Loss_t - Expenses_t) / Capital_t - \beta \cdot \Delta Reinsurance_t$$其中 $\alpha$ 是资本转化效率系数通常0.8~0.95$\beta$ 是再保险成本敏感度0.1~0.3。该模型比静态比率更有说服力因为它揭示了“一次大灾后需要几个季度恢复”。def simulate_car_evolution(initial_car2.1, n_quarters24, premium_growth0.015, loss_volatility0.25): 模拟资本充足率24个季度演化 initial_car: 初始CAR值如2.1表示210% premium_growth: 季度保费增长率小数如0.0151.5% loss_volatility: 损失波动率控制GPD抽样离散度 car_history [initial_car] capital_base 1000 # 初始资本基数百万美元 for t in range(n_quarters): # 保费收入按增长率递增叠加±5%随机扰动 premium 100 * (1 premium_growth)**t * (1 np.random.uniform(-0.05, 0.05)) # 损失支出从GPD抽样但限制单季损失不超过资本的30% loss_sample genpareto.rvs(shape, loc0, scalescale, size1)[0] u_threshold loss min(loss_sample, 0.3 * capital_base) # 防止单季清零 # 运营费用保费的22% expenses 0.22 * premium # 再保险变动当CAR1.8时自动增加再保分出成本上升15% reinsurance_delta 0.0 if car_history[-1] 1.8 else 0.15 # 更新CAR公式中α0.92, β0.25 delta_car 0.92 * (premium - loss - expenses) / capital_base - 0.25 * reinsurance_delta new_car car_history[-1] delta_car car_history.append(max(new_car, 0.5)) # CAR不低于50%否则视为技术性违约 return np.array(car_history) # 运行100次蒙特卡洛模拟 car_simulations np.array([simulate_car_evolution() for _ in range(100)]) car_mean np.mean(car_simulations, axis0) car_5pct np.percentile(car_simulations, 5, axis0) # 5%分位数监管关注底线 plt.figure(figsize(10, 5)) plt.plot(car_mean, b-, label平均CAR) plt.fill_between(range(25), car_5pct, alpha0.3, colorred, label5%分位数) plt.axhline(y1.5, colork, linestyle--, label监管阈值150%) plt.xlabel(Quarter) plt.ylabel(Capital Adequacy Ratio) plt.legend() plt.title(Capital Adequacy Ratio Evolution (100 Simulations)) plt.grid(True) plt.show()3. 将可持续性转化为可执行的费率优化策略建模的终点不是画一条CAR曲线而是回答“下季度对沿海地区砖混结构住宅基准费率应上调多少才能使CAR_5pct稳定在1.5以上” 这要求将统计模型与精算定价逻辑打通。我们采用梯度提升树XGBoost作为连接层输入是风险因子地理指数、建筑类型、免赔额、市场因子同业费率、通胀率输出是费率调整系数。关键创新在于损失函数不是常规的MSE而是CAR约束下的定制损失。3.1 构造带监管约束的XGBoost训练目标传统XGBoost回归会最小化预测误差但保险定价需满足调整后费率必须使模拟CAR ≥ 1.5。我们设计两阶段目标第一阶段用历史数据训练基础费率模型第二阶段在基础模型上叠加一个“约束校正项”该修正项由另一个XGBoost学习其标签是若当前费率导致CAR 1.5则标签为正需提价否则为0。import xgboost as xgb from sklearn.model_selection import train_test_split # 假设df_pricing为定价建模数据集含features列和target列历史费率 X df_pricing[feature_columns] y df_pricing[rate_per_thousand] # 第一阶段基础费率模型 X_train, X_test, y_train, y_test train_test_split(X, y, test_size0.2, random_state42) base_model xgb.XGBRegressor(n_estimators100, max_depth4, learning_rate0.1) base_model.fit(X_train, y_train) # 第二阶段约束校正模型 —— 标签为二元1需提价以满足CAR≥1.50无需 # 关键标签生成需调用2.3节的CAR模拟函数 def generate_constraint_label(row): # 用row中的风险因子模拟CAR若5%分位CAR1.5则返回1 simulated_car simulate_car_evolution( initial_car2.0, premium_growth0.015, loss_volatility0.25 ) return 1 if np.percentile(simulated_car, 5) 1.5 else 0 # 为训练集每行生成标签实际中需批量计算此处简化 y_constraint [generate_constraint_label(row) for _, row in X_train.iterrows()] constraint_model xgb.XGBClassifier(n_estimators50, max_depth3, learning_rate0.15) constraint_model.fit(X_train, y_constraint) # 预测时基础费率 × (1 0.15 × constraint_prediction) y_pred_base base_model.predict(X_test) y_pred_constraint constraint_model.predict(X_test) final_rate y_pred_base * (1 0.15 * y_pred_constraint)3.2 可视化费率调整的地理热力图与可持续性收益MCM 论文要求空间可视化。我们用 GeoPandas 绘制 ZIP Code 级别费率调整幅度并叠加 CAR 改善效果import geopandas as gpd import contextily as ctx # 加载美国ZIP Code地理边界需提前下载shp文件 gdf_zip gpd.read_file(tl_2023_us_zcta520.shp) gdf_zip gdf_zip.to_crs(epsg3857) # 转为Web Mercator投影 # 合并费率调整数据假设df_rates含zip_code和rate_change_pct列 gdf_merged gdf_zip.merge(df_rates, left_onZCTA5CE20, right_onzip_code, howleft) fig, ax plt.subplots(1, 2, figsize(18, 8)) # 左图费率调整热力图 gdf_merged.plot(columnrate_change_pct, cmapRdYlBu_r, legendTrue, axax[0], missing_kwds{color: lightgrey}) ax[0].set_title(Rate Adjustment by ZIP Code (%)) ctx.add_basemap(ax[0], sourcectx.providers.OpenStreetMap.Mapnik) # 右图CAR改善对比调整前vs调整后 car_before car_simulations[0] # 任选一次模拟 car_after simulate_car_evolution( initial_car2.0, premium_growth0.015 0.003, # 假设整体费率提升0.3% loss_volatility0.25 ) ax[1].plot(car_before, r--, labelCAR before adjustment) ax[1].plot(car_after, g-, labelCAR after adjustment) ax[1].axhline(y1.5, colork, linestyle:, labelRegulatory threshold) ax[1].set_xlabel(Quarter) ax[1].set_ylabel(Capital Adequacy Ratio) ax[1].legend() ax[1].grid(True) ax[1].set_title(CAR Evolution: Impact of Rate Adjustment) plt.tight_layout() plt.show()提示contextily的 basemap 需联网加载若离线环境可替换为gdf_zip.boundary.plot(axax[0])仅绘制行政边界。4. 在MCM限时场景下快速验证模型可靠性的3个硬核技巧竞赛中没有时间做全量交叉验证。我们用三个可在30分钟内完成的实证检验直接回应评审最可能质疑的点模型是否过拟合参数是否合理结论是否稳健4.1 用“反事实删除法”检验关键参数敏感性评审常问“如果GPD的shape参数变化±0.1VaR会偏移多少” 我们不重新拟合100次而是用解析法快速计算# GPD的VaR对shape的导数解析解 def var_sensitivity_to_shape(shape, scale, u, p0.995): 计算VaR_p对shape的偏导数 term (1 - p)**(-shape) d_term_d_shape -term * np.log(1 - p) return (scale / shape**2) * (term - 1) - (scale / shape) * d_term_d_shape sens var_sensitivity_to_shape(shape, scale, u_threshold) print(fVaR_99.5%对shape的敏感度: {sens:.2f} (每单位shape变化导致VaR变化)) # 若sens 50说明VaR对shape极度敏感需在论文中讨论shape的贝叶斯估计4.2 用Bootstrap重采样验证CAR阈值的统计显著性判断“CAR_5pct是否真的≥1.5”不能只看点估计。我们用Bootstrap计算其95%置信区间from sklearn.utils import resample def bootstrap_car_ci(simulations, alpha0.05): 对CAR序列的5%分位数计算置信区间 n_sim simulations.shape[0] n_boot 1000 boot_5pct [] for _ in range(n_boot): # 从100次模拟中重采样100次允许重复 indices np.random.randint(0, n_sim, n_sim) boot_sim simulations[indices] boot_5pct.append(np.percentile(boot_sim[:, -1], 5)) # 取最终季度的5%分位 lower np.percentile(boot_5pct, alpha/2 * 100) upper np.percentile(boot_5pct, (1 - alpha/2) * 100) return lower, upper lower_ci, upper_ci bootstrap_car_ci(car_simulations) print(fCAR_5pct at quarter 24: [{lower_ci:.3f}, {upper_ci:.3f}]) # 若区间完全在1.5上方可写“在95%置信水平下可持续性成立”4.3 构建“压力测试仪表盘”一键生成关键图表把上述所有检验封装成函数输入新数据即可输出PDF报告def generate_stress_dashboard(new_data_path): 输入新赔案CSV输出含4张图的PDF df_new pd.read_csv(new_data_path) # ... 执行2.1~2.3节全部流程 ... # 生成四图1. Mean Excess Plot 2. CAR模拟分布 3. 费率热力图 4. 敏感性矩阵 fig, axes plt.subplots(2, 2, figsize(16, 12)) # ... 绘图代码 ... plt.savefig(stress_test_report.pdf, bbox_inchestight) print(压力测试报告已生成: stress_test_report.pdf) # 竞赛中只需一行命令 # generate_stress_dashboard(data/2024Q1_claims.csv)注意bbox_inchestight可自动裁剪图表白边避免PDF中出现大片空白这是MCM排版的隐形加分项。本文还有配套的精品资源点击获取