MATLAB：数字期权定价、蒙特卡洛与显式积分公式？

Question

我在使用 MATLAB 时遇到以下问题：

令 Z 服从对数正态分布，使得 ln Z 具有均值 m 和方差 w。设 eta 为负数，c 为正常数。

我正在尝试计算期望值（令 I(Z<=c) denote the indicator function of the set (Z<=c))

E[Z^(eta+1) I(Z<=c)] = (1/sqrt(w)) integral_0^c x^(eta) phi((ln x - m)/sqrt(w)) dx,

其中 phi() 表示标准正态随机变量的概率分布函数。

我做的第一件事是模拟 10.000 次 Z 试验，将值 >c 的向量条目设置为 0，计算 (eta+1) 次方，然后计算平均值。这应该给我预期值的 MC 估计。

ST = random('Lognormal', m,w_sq,10000,1);
hlp = zeros(10000,1);
hlp(ST<=2) = ST(ST<=2);
hlp(hlp>0) = hlp(hlp>0).^(eta1+1); % 0^(eta1+1) gives infinity
mean(hlp)

对于积分，我使用了以下代码

tmpp = integral(@(x) x.^(eta1) .* normpdf((log(x)-m)/sqrt(w_sq),0,c);
tmpp / sqrt(w_sq(1))

不幸的是，这些过程会导致完全不同的结果，尽管从数学上讲它们应该是相同的。

这整个事情是一个更大代码的一部分，对我来说使用完整版本会方便得多。最初我尝试使用 MC 模拟进行双重检查，然后发现一定有问题......

有人可以帮忙吗？

Answer 1

对于第一段代码，

我想你之所以能得到意想不到的结果，是因为当你 对

hlp

进行计算时，您会尽量避免使用 0 值，如

0^eta

会爆炸 - 这不是想要的结果，所以你只需放弃它即可。 但在最后一步中，

mean(hlp)

将获取数组

hlp

中的所有值， 包括那些0。试试这个：

mean(hlp(hlp>0))

.

我的结果大约是 2.x，10,000 点模拟，2.3x ~ 2.4 倍，1,000,000 点。

我错了。你的问题是你使用的点数太少了。试试10,000,000积分，你会满意的:)

对于第二个问题，我无法理解定义变量的方式。 （我没有足够的声誉来添加评论，所以我把它放在这里。）

w_sq

中的“sq”是否意味着

w

的平方根？因为根据

random

的文档，参数应该是“sigma”，即标准差。将 SD 定义为方差的平方根是很自然的，我猜它是

w

。那么你在第一部分就做得对了。

如果是这样，在你的第二段代码中，为什么要把

sqrt()

放在

w_sq

上？你的意思是求方差的四次方根吗？从你对期望的定义来看，我认为这是不正确的。请看一下。

另一方面，

w_sq

是单个数字还是数组？

• 如果是数字 -

tmpp / sqrt(w_sq(1))

应该是

tmpp / sqrt(w_sq)

，尽管事实上它们没有什么区别。

• 如果它是一个数组 -

您可能希望将所有代码放入

for

循环中。循环会迭代

w_sq

，每次它都会挑选出数组中的一个元素（将变量命名为

w_sq_elem

），并让其余代码像单个数字一样执行操作。

总之，

(log(x)-m)/sqrt(w_sq)

和

tmpp / sqrt(w_sq(1))

讲述的是关于

w_sq

的不同信息。第一个假设它是一个数字，因此除法可以简单地是

/

，而不是

./

。第二个表明它是一个数组，因此您正在选择它的第一个元素。但数组在这里没有意义，因为根据我的理解，你并不是将

x

中的 10,000 个点除以 10,000 个不同的方差。

m = 3;
w_sq = 2;
eta1 = -2;
ST = random('Lognormal',m,w_sq,10000000,1);
hlp = zeros(10000000,1);
hlp(ST<=2) = ST(ST<=2);
hlp(hlp>0) = hlp(hlp>0).^(eta1+1); % 0^(eta1+1) gives infinity
mean(hlp)

tmpp = integral(@(x) x.^(eta1) .* normpdf((log(x)-m)/w_sq),0,2) ;
tmpp / w_sq

>> untitled

ans =

    0.2944


ans =

    0.2948

>>