尧图网站设计 尧图网站设计YAOTU DESIGN
ARTICLE DETAIL

资讯详情

深耕网站设计与一线实操的经验洞察。

MATLAB中的马尔可夫链与风险博弈建模:从状态矩阵到吸收概率

MATLAB中的马尔可夫链与风险博弈建模:从状态矩阵到吸收概率 做风险决策的时候最常被问的一句话是“还要不要再投入一次”表面看这是个凭胆量拍板的时刻但拆开细看你会发现它其实是一个典型的随机过程——每一次投入的收益、亏损、持平都像一个带概率的“下一步”而下一步会走到哪里只取决于你当前手里还剩多少筹码。这个“下一步只依赖当前状态”的性质正是马尔可夫链的核心。我最近在MATLAB里把这类“风险博弈”完整建模并求解了一遍从状态转移矩阵到吸收概率再到仿真验证整个过程走下来收获很大这篇就把它记录下来。无论你是学运筹学、做风险管理还是单纯对随机过程建模感兴趣这套“从规则到矩阵再到决策”的分析链路都值得抄走。1. 风险博弈里的马尔可夫性藏在哪几个关键假设里1.1 一场风险博弈怎么变成一个状态序列很多人听到“马尔可夫链”第一反应是数学公式、矩阵乘法、特征值分解觉得离实际问题很远。但恰恰相反风险博弈是最适合用马尔可夫链去描述的日常场景之一。你可以把一个参与者当前的资金、积分、资源量、剩余寿命统统抽象成一个“状态”。每一次决策或者每一轮随机事件发生后状态从一个数值跳到另一个数值这种跳变只跟当前状态以及当前规则有关跟你过去是怎么走到这里的没有直接关系。这个“无记忆性”听起来好像很苛刻但现实中大量风险模型恰好满足。举个最常见的例子我帮朋友整理资金管理模型时他把账户余额作为唯一关心的变量。单次博弈的三种结果——盈利一个单位、亏损一个单位、持平——每轮独立同分布下一轮余额只取决于当前余额是多少至于上一轮是盈利还是亏损并不会改变下一轮的概率分布。这就是马尔可夫性。从数学角度表述就是设随机变量序列为 (X_0, X_1, X_2, \dots)每个取值代表一个状态当对所有 (k) 都有[ P(X_{k1}j \mid X_ki) P(X_{k1}j \mid X_ki, X_{k-1}i_{k-1}, \dots, X_0i_0) ]即向前转移的概率只依赖当前状态 (i) 而不依赖更早的历史时这个序列就可以用状态转移矩阵来描述。对于风险博弈来说(X_k) 最自然的选择就是“第 (k) 轮后的剩余资源量”。很多人第一次接触马尔可夫链是通过“醉汉随机游走”模型一个醉汉在街上随机迈步走得再远下一步只取决于他现在站在哪里。风险博弈里的资金变动本质上就是带边界的随机游走——游走者有可能会踩到“0”这个边界资源耗尽也有可能走到某个目标边界达成目标收益。所以酒精模型和风险管理模型在数学上是同一类东西。1.2 经典Gamblers Ruin模型从一个历史问题到一般风险分析“Gamblers Ruin”这个名字听起来有点赌场味道但它在学术和应用层面是个非常经典的马尔可夫链问题。模型背景很朴素参与者有初始资金 (i)每轮以概率 (p) 赢一个单位、以概率 (q1-p) 输一个单位当资金涨到目标值 (N) 时博弈结束或者跌到 0 时宣告破产。问题是从 (i) 出发最终到达 (N) 而不是 0 的概率是多少平均要经历多少轮才能结束如果只把它理解成“一个赌徒会不会输光”那确实太窄了。这个模型的本质是一套“带吸收边界的随机游走求解工具”在风险管理里的应用极其广泛保险公司现金流在“偿付能力充足”和“破产”之间波动量化策略的账户权益在“到达止盈线”和“触发止损线”之间上下跳动甚至设备运维里“正常运行”和“故障停机”两个吸收态之间的转换也可以用同一套模型来描述。关键不在于名字而在于你能不能把实际问题抽象成“几个状态 一套转移概率”。用MATLAB来做这件事最大的优势就是矩阵运算和仿真代码都能在同一个环境里完成。你不需要手动去解复杂的差分方程只要构建出转移矩阵剩下交给线性代数库就好。下面我从建模开始一直讲到最后怎么把概率结果翻译成实际决策。2. 状态空间定义与转移矩阵构建建模是最容易出错的一步2.1 状态空间的闭合性0和N为什么必须是吸收态无论是解析求解还是MATLAB仿真第一步永远是定义状态空间。在Gamblers Ruin场景中状态集可以设为 ({0, 1, 2, \dots, N})其中 (0) 和 (N) 是吸收态。所谓吸收态就是一旦进入就再也出不来的状态资金变成 0游戏结束无法再继续资金达到目标 (N)同样结束博弈。这两个状态在转移矩阵里体现为“自己转移到自己的概率为 1”。这个“闭合性”特别重要是新手最容易忽略的点。很多人会下意识地认为“资金达到0之后可能还会有人借钱继续”或者“达到N之后还能再继续玩”如果你的实际问题确实是那样那状态空间就要重新定义。但在一类“到达边界即停止”的风险博弈里边界必须是吸收态否则马尔可夫链无法收敛求解“最终到达某个边界的概率”就无从谈起。定义好吸收态之后还需要验证状态空间的完备性从任意一个中间状态出发经过一步转移到达的所有状态必须仍然落在 ({0,1,\dots,N}) 里面。如果规则里存在“一次性赢3个单位”或者“一次性输5个单位”的情况那么状态空间就不能只按单位资金粒度来建可以考虑把每个可能到达的余额都列出来或者调整状态粒度。边界条件写错了后面矩阵再漂亮也算不出正确结果。2.2 从规则到矩阵转移概率矩阵的完整推导假设每轮只有三种结果以概率 (p) 赢一单位、以概率 (q) 输一单位、以概率 (r) 持平并且满足 (pqr1)。那么对于任意中间状态 (i \in {1,2,\dots,N-1})转移规则就是[ P_{i,i1} p,\quad P_{i,i-1} q,\quad P_{i,i} r ]对于边界有[ P_{0,0} 1,\quad P_{N,N} 1 ]把 (N1) 个状态按 (0,1,\dots,N) 的顺序排列转移矩阵就是一个 ((N1)\times(N1)) 的三对角矩阵。对角线上的平局概率 (r) 会影响矩阵性质吗会但不会改变吸收概率的最终结果。因为平局只是不改变资金数相当于在原地多等了一轮真正决定你是不是能到边界的还是 (p) 和 (q) 的比例关系。这里有一个很实用的建模建议与其直接用完整矩阵不如先把状态分成“暂态”和“吸收态”两个部分。暂态集合是 (T {1,\dots,N-1})吸收态集合是 (A {0, N})。把转移矩阵写成标准块结构[ P \begin{bmatrix} Q R \ 0 I \end{bmatrix} ]其中 (Q) 是暂态之间的转移概率矩阵(R) 是暂态转移到吸收态的概率矩阵(0) 是吸收态到暂态的零矩阵(I) 是吸收态内部的自转移矩阵。这种写法是后续求吸收概率和期望吸收时间的基础也是MATLAB里高效计算的起点。2.3 一个可直接运行的示例目标金额与止损线设定为了后面所有代码都能跑起来我先给一个具体的案例参数。假设参与者初始资金为 100每次按固定单位下注目标收益是涨到 150即 (N150)止损线是跌到 50 就强制离场。数学上为了统一可以把状态平移让止损线对应状态 0目标线对应状态 100。但为了方便理解我在代码里直接用 (N150)、下边界 (0)、初始状态 (100) 来做。状态集就是 ({0,1,\dots,150})。单轮概率先设为 (p0.48, q0.46, r0.06)也就是每轮有48%概率盈利一个单位46%概率亏损一个单位6%概率原地踏步。这个参数设置演示的是“单轮胜率略低于50%但平局缓冲了一部分风险”的场景在现实里比如某些保险合同定价、库存损耗模型中都有对应含义。用MATLAB构建转移矩阵的代码如下N 150; p 0.48; q 0.46; r 0.06; init 100; % 构建完整转移矩阵 P P zeros(N1, N1); P(1,1) 1; % 状态0破产吸收态 P(N1,N1) 1; % 状态N目标吸收态 for i 2:N P(i, i1) p; % 从状态 i-1 赢到 i P(i, i-1) q; % 从状态 i-1 输到 i-2 P(i, i) r; % 平局 end这里下标为什么要从2到N而不是1到N因为MATLAB数组从1开始状态0对应数组索引1状态1对应索引2以此类推。这个映射关系虽然琐碎但写代码的时候只要错一位整个矩阵就全乱了。建议在代码开头用注释把“状态值 ↔ 数组索引”的映射写清楚这也是我自己写这类模型时养成的习惯。3. MATLAB求解吸收概率与期望轮数解析解与蒙特卡洛双重验证3.1 用线性方程组的观点避免矩阵求逆状态转移矩阵构建完之后下一个问题就是怎么求“从初始状态出发最终到达目标吸收态的概率”。很多人第一反应是对 (P) 做矩阵连乘看极限分布。在状态数少的时候这确实可行但状态数一多连乘的效率很低而且数值上容易出问题。更规范的解法是直接从吸收马尔可夫链的标准理论出发求解线性方程组。设 (h_i) 表示从状态 (i) 出发最终到达吸收态 (N) 的概率那么对所有暂态 (i) 有[ h_i \sum_{j \in T} P_{i,j} h_j \sum_{k \in A} P_{i,k} h_k ]写成矩阵形式[ h_T Q h_T R h_A ]由于吸收态 (0) 的最终到达目标概率是0吸收态 (N) 的最终到达目标概率是1所以 (h_A [0, 1]^T)。整理得到[ (I - Q) h_T R \cdot \begin{bmatrix} 0 \ 1 \end{bmatrix} ]这个线性方程组只需要解一次MATLAB直接用反斜杠运算符\就能搞定不需要手写任何迭代。反斜杠运算符在MATLAB里是对线性方程组求解的高度优化实现对中小规模问题稳定且高效。期望吸收轮数的方程同理。设 (e_i) 表示从状态 (i) 出发到最终吸收所需的期望轮数则有[ e_T \mathbf{1}_T Q e_T ]其中 (\mathbf{1}_T) 是全1列向量意思是每一轮先消耗1步然后再加上后续的期望步数。整理后[ (I - Q) e_T \mathbf{1}_T ]3.2 解析解实现p≠q和pq两种场景基于块结构MATLAB里可以这样实现暂态索引的提取transientIdx 2:N; % 对应状态 1 到 N-1 absorbingIdx [1, N1]; % 对应状态 0 和 N Q P(transientIdx, transientIdx); R P(transientIdx, absorbingIdx); % 目标吸收态 N 在 absorbingIdx 中是第2列 b R(:, 2); I eye(size(Q)); h (I - Q) \ b; % h(transient_state_index) % 初始状态 init100 在 transientIdx 中的位置是 100 initTransIdx init; probReachTarget h(initTransIdx); % 期望吸收轮数 e (I - Q) \ ones(size(Q,1), 1); expectedRounds e(initTransIdx);运行后你会得到一个具体的数字比如当 (p0.48, q0.46) 时初始资金100到目标150的概率大约只有1.8%期望轮数则可能高达几千轮。这个结果非常反直觉——单轮胜率只差2个百分点最终实现目标的概率却发生崩塌式下降。这就是马尔可夫链模型的价值它把单轮概率的微小偏差通过长期链条放大成了最终结果的数量级差异。如果你的模型里 (pq0.5, r0)解析解有一个经典闭合公式从状态 (i) 出发到状态 (N) 的概率是 (i/N)期望轮数是 (i(N-i))。这个公式非常适合用来验证代码逻辑。我在调试时通常先把参数设成这种简单场景确认MATLAB输出等于理论值后再换成真实参数。这个习惯帮我省下了大量排查矩阵索引错误的时间。3.3 蒙特卡洛仿真验证解析解的“照妖镜”解析解算出来之后我强烈建议再写一个蒙特卡洛仿真去交叉验证。原因很简单解析解依赖矩阵方程的正确性一旦矩阵某一行写错结果照样是错的但你可能看不出来。仿真代码的逻辑更贴近原始规则如果两边的结果能对上说明建模到求解的整个链路都没问题。下面是蒙特卡洛仿真代码trials 20000; countWin 0; totalSteps 0; for t 1:trials state init; step 0; while state 0 state N x rand(); if x p state state 1; elseif x p q state state - 1; end % 剩下概率 r 为平局state 不变 step step 1; if step 1e6 break; end end if state N countWin countWin 1; end totalSteps totalSteps step; end fprintf(蒙特卡洛达标概率: %.4f\n, countWin / trials); fprintf(蒙特卡洛平均轮数: %.2f\n, totalSteps / trials); fprintf(解析解达标概率: %.4f\n, probReachTarget); fprintf(解析解平均轮数: %.2f\n, expectedRounds);我在实际运行中20000次试验得到的达标概率大概在1.7%到2.0%之间浮动解析解在1.8%附近非常吻合。平均轮数的方差比较大可能需要更多试验次数才能稳定到解析解的水平但趋势是对的。有一点要提醒蒙特卡洛只是验证工具不是求解主工具。在状态数很大、目标概率很小比如万分之一的情况下仿真需要海量试验次数才能看到一次达标事件效率极低。而解析解法无论概率多小只要矩阵能解结果都是精确的。所以我个人的实践原则是小规模问题用仿真验证正式计算用解析解。4. 把数值结果翻译成决策动作收益目标与止损线的权衡4.1 解读吸收概率你的目标不是“赢一次”是“到边界”模型跑完之后最大的挑战不是拿到概率而是理解这个概率对决策意味着什么。很多人的直觉是“单轮赢的概率是48%那长期下来应该接近48%的胜率吧”这个想法在无边界随机游走里可能有一定道理但一旦加了止损线和目标线事情就完全不同了。当前参数下从100到150的概率只有1.8%左右而到50触发止损的概率是98.2%。为什么差距这么大原因是这个随机游走带有明显的负漂移——每个单位资金的期望变化是 (p \times (1) q \times (-1) r \times 0 0.48 - 0.46 0.02)。等等这样算期望变化其实是正的0.02为什么实现目标的概率还那么低这个计算需要小心虽然单轮期望值是正的但目标距离和止损距离在这里正好相等都是50个单位而单轮输的概率更大随机波动会导致链条末端更密集地聚集在较近的“劣势”区域。更直观的理解方式是这样的从100出发上升需要连续积累较多正向净胜但一旦产生向下波动离止损线50更近。因为 (q p)系统有明显的下降趋势哪怕每轮优势看起来只差2个百分点在50步需要净赢这个量级的要求下这种劣势会被指数放大。所以结论是在目标距离和止损距离相当时任何单轮负优势都会导致达标概率迅速塌缩到接近零。单轮胜率 p单轮败率 q达标概率解析解期望轮数0.500.500.666750000.490.490.1108约47000.480.460.0183约45000.470.530.0025约4300这张表就是“把数值结果翻译成决策”的核心单轮胜率从0.50降到0.48实现目标收益的概率从三分之二掉到不到2%。它传达的决策含义很清晰——如果你的策略长期单轮胜率没法稳定超过50%那么设置“距离相同的目标线和止损线”是一个非常危险的架构。4.2 参数敏感性胜率变动5个点结果差多少用MATLAB做参数扫描非常方便只需要把上面的求解逻辑包进一个循环遍历不同的 (p) 值即可。比如对 (p \in {0.45, 0.46, 0.47, 0.48, 0.49, 0.50}) 逐一求解并画出达标概率的变化曲线。实际跑下来你会发现这段曲线几乎是“指数型”的从0.5附近的较高概率一路下滑到0.45时的几乎为0。为什么会出现这种非线性因为在吸收概率公式里真正起作用的是比值 (q/p) 的幂次。初始资金100、目标150这个设置意味着链条要经历足够多的随机波动而 (q/p) 只要偏离1一点点经过几十次方之后就会被放大成巨大的倍数直接压垮向上到达目标的概率。这种现象叫“尾部风险放大”用马尔可夫链模型看是最直观的。如果你做的是实际策略评估这个敏感性分析应该做的不是自己拍脑袋而是把所有关键参数都纳入分析。比如赢利幅度与亏损幅度不对称赢一单位、输两单位这时单纯用胜率分析就不够了需要把状态转移的步长改成非对称形式但求解框架完全不变。MATLAB的循环 矩阵求解在这种场景下效率非常高扫一组参数也就是几十毫秒。4.3 用同一框架比较不同博弈规则同一个马尔可夫链框架最大的价值是可以横向比较不同的博弈规则。比如你设计了两种方案方案A是目标150、止损50、单轮胜率0.48方案B是目标120、止损80、单轮胜率0.48。把 (N) 改一下重新求解你会发现方案B的达标概率可能显著高于方案A。这说明在胜率不变的情况下调整目标线和止损线的相对位置是改善风险收益结构的有效手段。再比如方案C同样是目标150、止损50但规则改成“盈利时下一次盈利概率提高到0.55亏损后下一次盈利概率降到0.40”这个时候严格来说已经不是标准一阶马尔可夫链了因为转移概率依赖上一轮的结果。但你可以把状态空间扩成“盈亏方向×资金水平”的联合状态模型依然可以用相同框架求解。状态空间变大带来的代价是矩阵维度增加但MATLAB的稀疏矩阵能轻松处理几万维的问题。如果你做的项目需要出决策报告我建议把每个候选方案都输出一张表格达标概率、止损概率、期望轮数、概率分布的另一侧分位数。这些指标全部可以靠上面的代码扩展得到不需要额外引入复杂的计算工具。5. 这套模型的边界、常见坑与扩展思路5.1 马尔可夫性假设何时失效记忆效应与动态策略我在前面不断强调“无记忆性”但这个假设不是什么时候都成立。如果决策者的策略本身依赖历史信息——比如“连续亏损两次后降低仓位”“连续盈利三次后提前止盈”——那么当前余额就不足以描述未来演化的概率必须把历史信息也纳入状态。一个常见的做法是构造扩展状态把“上次是否盈利”加进状态编码形成类似 ((资金水平, 上一轮结果)) 的组合状态。这样做之后转移矩阵依然满足马尔可夫性只是维度翻倍求解思路不变。另一种失效场景是外部环境随时间变化。如果单轮胜率不是常数而会随着宏观条件变化那就成了非齐次马尔可夫链。非齐次模型更难处理但实务中可以先分段处理把时间切成多个区间在每个区间内假定转移概率不变然后分段相乘或者分段求解。这属于工程近似不是严格数学推导但应对实际问题已经足够。5.2 MATLAB实现中的数值与实践问题用MATLAB做这类模型有几个坑我必须提一下。第一矩阵维度过大时不要用全矩阵。当 (N10000) 时((N1)^2) 大约是1亿个元素double类型占800MB内存很容易把MATLAB卡死。必须改用sparse稀疏矩阵来构建转移矩阵。好消息是转移矩阵天然是稀疏的三对角结构在稀疏存储下只占 (O(N)) 内存。N 10000; p 0.48; q 0.46; r 0.06; P sparse(N1, N1); P(1,1) 1; P(N1,N1) 1; P(2:N, 3:N1) p; P(2:N, 1:N-1) q; P(2:N, 2:N) r;这段代码用了稀疏矩阵的向量化赋值避免for循环在大 (N) 时构造速度会快非常多。注意索引边界P(2:N, 3:N1)p表示从状态 (i) 到 (i1) 的转移这里行索引2到N对应状态1到N-1列索引3到N1对应状态2到N。第二(I-Q)\b在 (N) 很大的时候也可能变慢。这时可以考虑使用迭代求解器比如bicgstab、gmres或者直接利用三对角结构写专用的Thomas算法。不过对于几千到几万的状态规模MATLAB自带的\已经很快我用下来基本没有卡顿。第三蒙特卡洛仿真要注意随机数种子。调代码阶段最好设置rng(0)保证每次运行结果一致方便比对。等正式跑研究结果时再取消固定种子避免结论过度依赖某一次随机序列。5.3 扩展方向高阶马尔可夫链与运筹学工具链这个模型做完之后我最大的体会是马尔可夫链不是孤立存在的数学玩具它是整个运筹学决策框架里的基础零件。在《运筹学基础及其MATLAB应用》这类教材里你会看到马尔可夫链跟决策树、动态规划、排队论反复交叉出现。比如说你可以在马尔可夫链算出的达标概率基础上再叠加一个“是否值得投入”的期望效用函数也可以把吸收概率作为动态规划的状态转移收益进一步设计最优的止盈止损策略。我自己在项目里往后扩展过几个方向一是把平局概率 (r) 从常数改成跟当前状态相关比如接近止损线时平局概率升高这种“近边界惰性”在真实系统中很常见二是引入时间成本给期望轮数加上贴现因子让“更早到达”和“更高概率到达”之间形成可比较的指标三是把单变量状态扩展成二维状态比如“资金 × 时间裕度”分析有期限限制下的最优博弈策略。如果你也打算把这套模型写进论文或报告我强烈建议每次跑完解析解之后都顺手跑一组蒙特卡洛仿真。别嫌多此一举这个习惯救过我很多次——有时候矩阵里的索引偏移一位解析解仍然能输出一个看似合理的数字只有仿真能帮你抓住这种隐藏极深的错误。模型能不能从纸面走向落地靠的不是某个高深算法而是这些老老实实的验证习惯。
返回列表