MATLAB：求解一系列方程以系统地调整概率估计答案

【问题标题】：MATLAB: Solving series of equations to systematically adjust probability estimatesMATLAB：求解一系列方程以系统地调整概率估计
【发布时间】：2023-03-05 02:40:01
【问题描述】：

我正在尝试编写 MATLAB 代码，接收估计的概率 x 和 y 并将它们映射到一组调整后的概率 p 和 q。

通常输入概率x 和y 不会加起来为 1，因此代码的重点是调整概率，使它们加起来为 1。

这是一个例子，x + y = 1.03

syms n p q x y
eqn1 = x == 0.7813;
eqn2 = y == 0.2487;
eqn3 = p^n == x;
eqn4 = q^n == y;
eqn5 = log(x)./log(p) == log(y)./log(q);
eqn6 = x^(1/n) == p;
eqn7 = y^(1/n) == q;
eqn8 = p + q == 1;

sol = solve([eqn1, eqn2, eqn3, eqn4, eqn5, eqn6, eqn7, eqn8], [n p q x y])

n 的值同时应用于 x 和 y，它故意比高概率估计更能“压低”低概率估计。

我知道在上面的示例中，答案应该是n≈0.946139, p≈0.7704, q≈0.2297，但无法编写任何产生这种输出的代码。

上面的代码产生一个空结构。我也尝试过vpasolve，但也没有用。我不确定我的解决方案是否失败是因为我错误指定了问题，还是因为我在 MATLAB 中使用了错误的函数类型。

【问题讨论】：

当solve返回一个没有警告的空结构时意味着没有解决方案，如果有警告说明解决方案没有找到解决方案但存在解决方案。你收到警告了吗？

标签： matlab

【解决方案1】：

这里的方程式有很多冗余。

eqn6 和 eqn7 分别暗示 eqn3 和 eqn4。此外，eqn3 和eqn4 暗示eqn5。我们也不需要将x 和y 定义为常量的eqn1 和eqn2，因为使用直接替换可以更有效地处理它们。我猜想这些冗余是 MATLAB 感到困惑的原因，但我不确定。

此外，为此使用符号工具箱可能有点矫枉过正，因为可能没有封闭形式的解决方案，但如果你必须这样做，那么这里有一个解决方案。

x = 0.7813;
y = 0.2487;
syms n p q
eqn6 = x^(1/n) == p;
eqn7 = y^(1/n) == q;
eqn8 = p + q == 1;
sol = solve([eqn6, eqn7, eqn8], [p q n]);

p = sol.p
q = sol.q
n = sol.n

结果：

Warning: Cannot solve symbolically. Returning a numeric approximation instead. 
> In solve (line 304)
  In q51941792_2 (line 7) 
p = 0.77034051295869623476715256144042
q = 0.22965948704130376523284743855958     
n = 0.94585921162748018347665683036153

请注意，求解器找不到封闭形式的解，因此它默认采用数值优化方法。

如果您有优化工具箱，那么您可以改用fsolve，这样更合适（而且速度更快）。要使用fsolve，您需要定义一个函数，该函数接受一个输入向量并返回一个输出向量。 fsolve 的目标是找到导致所有零输出的输入向量。我们重新排列eqn6、eqn7 和eqn8 来形成这个函数。您可以阅读文档以获取更多信息。

x = 0.7813;
y = 0.2487;

% define optimization function assumes ary = [p, q, n]
f = @(ary) [x^(1/ary(3)) - ary(1); ...   % eqn6
            y^(1/ary(3)) - ary(2); ...   % eqn7
            ary(1) + ary(2) - 1];        % eqn8

% initial estimate
p0 = x;
q0 = y;
n0 = 1;

% solve
fsolve_opts = optimoptions('fsolve', 'OptimalityTolerance', 0, 'FunctionTolerance', 1e-12);
sol = fsolve(f, [p0; q0; n0], fsolve_opts);

p = sol(1)
q = sol(2)
n = sol(3)

结果：

p = 0.770340512958696
q = 0.229659487041304
n = 0.945859211627480

存在解决方案

有些情况是不存在解决方案的。例如，如果x > 1 和y > 0 & y < 1 则p 和q 都是正数，但其中一个大于1（取决于n 的值），但这与eqn8 相矛盾，因此不存在解决方案。如果您想使用这种方法，您需要保证存在解决方案。我相信一个充分条件是x 和y 都在 0 和 1 之间，虽然我还没有证明。

【讨论】：