主要内容

使用GPU ARRAYFUN进行蒙特卡罗模拟

这个例子展示了如何使用蒙特卡罗方法在GPU上计算金融期权的价格。本文以三种简单的奇异期权为例,但更复杂的期权也可以用类似的方式定价。

这个例子是一个函数,因此帮助器可以嵌套在其中。

函数paralleldemo_gpu_optionpricing

本例使用的是长时间运行的内核,因此如果GPU上的内核执行超时,则无法运行。当选择的GPU也在驱动显示时,超时通常是有效的。

dev = gpuDevice ();如果dev.KernelExecutionTimeout错误(“pctexample: gpuoptionpricing: KernelTimeout”...如果GPU上的内核执行失败,这个例子将无法运行...“超时”。]);结束

股票价格演变

我们假设价格的演化符合对数正态分布,与无风险利率、股息收益率(如果有的话)和市场波动性相关。假设所有这些数量在选项的生命周期内都是固定的。由此得到价格的随机微分方程如下:

$${rm d}S = S \times
左\ [& # xA;(r d) {} \ rm d t + \σ\ε\ sqrt {{\ rm d} t} & # xA;\] $ $

在哪里年代美元是股票价格,r美元为无风险利率,$ d $是股票的年股息率,\σ美元是价格的波动和\ε美元表示高斯白噪声过程。假设美元(S + \δS) / S $为对数正态分布,可离散为:

$S_{t+1} = S_{t} × exp{
左\ [& # xA;\离开(r - d - \压裂{1}{2}\σ^ 2 \)\δt # xA;+根号{t}
正确\]& # xA;} $ $

举个例子,我们用100美元的股票每年1%的股息。中央政府利率假定为0.5%。我们检查了一个大致每天采样的两年时间窗口。市场波动性假定为每年20%。

上涨空间= 100;%股票价格100美元起。股息= 0.01;年股息率% 1%。riskFreeRate = 0.005;% 0.5%。timeToExpiry = 2;%选项的生存期,以年为单位。sampleRate = 1/250;%每年250个工作日。波动率= 0.20;% 20%波动。

我们重置随机数生成器以确保结果可重复。

种子= 1234;rng(种子);%重置CPU随机数生成器。gpurng(种子);重置GPU随机数生成器。

我们现在可以随着时间循环来模拟股票价格的路径:

价格=上涨空间;时间= 0;持有;time < timeToExpiry time = time + sampleRate;漂移= (riskFreeRate - dividend - volatility*volatility/2)*sampleRate;扰动=波动性*sqrt(sampleRate)*randn();价格=价格*exp(漂移+扰动);情节(时间、价格、“。”);结束;网格;包含(的时间(年));ylabel (的股票价格($));

在GPU上运行

为了在GPU上运行股票价格模拟,我们首先需要将模拟循环放入helper函数中:

函数finalStockPrice = simulateStockPrice(S,r,d,v,T,dT)t < t t = t + dT;dr = (r - d - v*v/2)*dT;pert = v*sqrt(dT)*randn();S = S*exp(dr + pert);结束finalStockPrice = S;结束

我们可以调用它数千次arrayfun.为了确保计算发生在GPU上,我们将输入价格设置为一个GPU矢量,每个模拟包含一个元素。为了准确地测量GPU上的计算时间,我们使用了gputimeit函数。

%创建输入数据。N = 1000000;startStockPrices =上涨空间*的(N, 1“gpuArray”);%运行模拟。finalStockPrices = arrayfun(@simulateStockPrice,...startStockPrices, riskFreeRate, dividend, volatility,...timeToExpiry sampleRate);meanFinalPrice =意味着(finalStockPrices);%使用gputimeit度量该功能在GPU上的执行时间。这需要我们将|arrayfun|调用存储在函数句柄中。functionToTime = @() arrayfun(@ simulatestockprice,)...startStockPrices, riskFreeRate, dividend, volatility,...timeToExpiry sampleRate);timeTaken = gputimeit (functionToTime);流('计算的平均价格$ 1.4f在%1.3f秒。\n'...meanFinalPrice timeTaken);clf;hist(最终股票价格,100);包含(的股票价格($)) ylabel (“频率”网格);
0.283秒内计算的平均价格为98.9563美元。

亚洲期权定价

作为一个例子,让我们使用一个基于股票在期权生命期内价格的算术平均数的欧式亚洲期权。我们可以通过在模拟过程中对价格进行累加来计算平均价格。对于看涨期权,如果平均价格高于执行价格,则执行期权,而支付是两者之间的差额:

函数optionPrice = asianCallOption(S,r,d,v,x,T,dT) T = 0;cumulativePrice = 0;t < t t = t + dT;dr = (r - d - v*v/2)*dT;pert = v*sqrt(dT)*randn();S = S*exp(dr + pert);累计价格=累计价格+ S;结束numSteps = (T / dT);meanPrice = cumulativePrice / numSteps;用今天的货币表示最后的价格。optionPrice = exp(-r*T) * max(0, meanPrice - x);结束

同样,我们使用GPU运行数千条模拟路径使用arrayfun.每个模拟路径给出了期权价格的独立估计,因此我们采用均值作为结果。

罢工= 95;%期权执行价格($)。optionPrices = arrayfun(@asianCallOption,...startStockPrices, riskFreeRate, dividend, volatility, strike,...timeToExpiry sampleRate);meanOptionPrice =意味着(optionPrices);%在GPU上测量执行时间并显示结果。@() arrayfun(@ asiancalloption,)...startStockPrices, riskFreeRate, dividend, volatility, strike,...timeToExpiry sampleRate);timeTaken = gputimeit (functionToTime);流('计算的平均价格$ 1.4f在%1.3f秒。\n'...meanOptionPrice timeTaken);
0.287秒内计算的平均价格为8.7210美元。

为回眸期权定价

对于这个例子,我们使用欧洲风格的回看期权,它的支付是期权生命周期内的最低股票价格和最终股票价格之间的差值。

函数optionPrice = euroLookbackCallOption(S,r,d,v,T,dT) T = 0;minPrice = S;t < t t = t + dT;dr = (r - d - v*v/2)*dT;pert = v*sqrt(dT)*randn();S = S*exp(dr + pert);如果S结束结束用今天的货币表示最后的价格。optionPrice = exp(-r*T) * max(0, S - minPrice);结束

请注意,在这种情况下,期权的执行价格是最低股票价格。因为最终股票价格总是大于或等于最小值,期权总是被执行,而不是真正的“可选”。

optionPrices = arrayfun(@ eurolookbackcallloption,...startStockPrices, riskFreeRate, dividend, volatility,...timeToExpiry sampleRate);meanOptionPrice =意味着(optionPrices);%在GPU上测量执行时间并显示结果。@() arrayfun(@ eurolookbackcallloption,)...startStockPrices, riskFreeRate, dividend, volatility,...timeToExpiry sampleRate);timeTaken = gputimeit (functionToTime);流('计算的平均价格$ 1.4f在%1.3f秒。\n'...meanOptionPrice timeTaken);
0.286秒内计算的平均价格为19.2711美元。

障碍期权的定价

最后一个例子使用了“向上并向外”障碍期权,如果股票价格达到障碍水平,该期权就失效。如果股票价格保持在障碍水平以下,那么最终的股票价格被用于正常的欧洲看涨期权计算。

函数optionPrice = upAndOutCallOption(S,r,d,v,x,b,T,dT) T = 0;(t < t) && (S < b) t = t + dT;dr = (r - d - v*v/2)*dT;pert = v*sqrt(dT)*randn();S = S*exp(dr + pert);结束如果S < b在障碍范围内,所以价格就像欧洲期权一样。optionPrice = exp(-r*T) * max(0, S - x);其他的%触及障碍,因此期权被撤回。optionPrice = 0;结束结束

请注意,我们现在必须提供期权的执行价格和期权失效的障碍价格:

罢工= 95;%期权执行价格($)。障碍= 150;%该选项的障碍价格($)。optionPrices = arrayfun(@ upandoutcallloption,...startStockPrices, riskFreeRate, dividend, volatility,...罢工,障碍,...timeToExpiry sampleRate);meanOptionPrice =意味着(optionPrices);%在GPU上测量执行时间并显示结果。@() arrayfun(@ upandoutcallloption,)...startStockPrices, riskFreeRate, dividend, volatility,...罢工,障碍,...timeToExpiry sampleRate);timeTaken = gputimeit (functionToTime);流('计算的平均价格$ 1.4f在%1.3f秒。\n'...meanOptionPrice timeTaken);
0.289秒内计算的平均价格为6.8166美元。
结束