線形計画法を使用した長期投資の運用最大化: 問題ベース

ライブスクリプトを開く

この例では、問題ベースのアプローチを使用して、固定の年数 T にわたって確定的なリターンがある投資の問題を解く方法を説明します。この問題は利用可能な投資に資金を割り当てて最終資産額を最大化するものです。ソルバーベースのアプローチについては、線形計画法を使用した長期投資の運用最大化: ソルバーベースを参照してください。

問題の定式化

N 個のゼロクーポン債に T 年の期間にわたって投資するための初期資金 Capital_0 があるとします。各債券の固定金利は年単位の複利で、満期時に元本と複利が支払われます。目的は T 年後に回収される総資産額を最大化することです。

投資をする時点で、単一の投資が総資本の特定の割合を超えないようにする制約を含めることができます。

この例では、まず小規模なケースで問題の設定を示し、次に一般的なケースに定式化します。

この問題は、線形計画問題としてモデル化できます。そこで、資産額を最適化するような、最適化問題アプローチを使用した問題を策定します。

基本的な例

まず、小規模な例から始めます。

初期投資額 Capital_0 は $1000 です。
期間 T は 5 年です。
債券数 N は 4 です。
未投資資金をモデル化するために、満期が 1 年、利率が 0% で毎年利用可能な 1 つのオプション B0 を設定します。
債券 1 は B1 で表され、1 年目に購入でき、満期は 4 年、利率は 2% です。
債券 2 は B2 で表され、5 年目に購入でき、満期は 1 年、利率は 4% です。
債券 3 は B3 で表され、2 年目に購入でき、満期は 4 年、利率は 6% です。
債券 4 は B4 で表され、2 年目に購入でき、満期は 3 年、利率は 6% です。

最初のオプション B0 を満期が 1 年で利率が 0% の 5 つの債券に分割することによって、利用可能な債券合計数が 9 であるとしてこの問題を等価的にモデル化できます。つまり、k=1..9 に対して、次のとおりとします。

ベクトル PurchaseYears のエントリ k は、債券 k を購入可能な年の始まりを表します。
ベクトル Maturity のエントリ k は、債券 k の満期 $m_{k}$ を表します。
ベクトル MaturityYears のエントリ k は、債券 k を売却可能な年の終わりを表します。
ベクトル InterestRates のエントリ k は、債券 k の利率 $ρ_{k}$ を表します。

各債券の購入可能回数と期間を水平のバーで表すことによって、この問題を視覚化します。

% Time period in years
T = 5;
% Number of bonds
N = 4;
% Initial amount of money
Capital_0 = 1000;
% Total number of buying oportunities
nPtotal = N+T;
% Purchase times
PurchaseYears = [1;2;3;4;5;1;5;2;2];
% Bond durations
Maturity = [1;1;1;1;1;4;1;4;3];
% Bond sale times
MaturityYears = PurchaseYears + Maturity - 1;
% Interest rates in percent
InterestRates = [0;0;0;0;0;2;4;6;6];
% Return after one year of interest
rt = 1 + InterestRates/100;

plotInvestments(N,PurchaseYears,Maturity,InterestRates)

決定変数

ベクトル x によって決定変数を表します。ここで、x(k) は、債券 k (k = 1,...,9) へのドル建て投資額です。満期時の投資 x(k) に対する支払いは次のようになります。

$x (k) (1 + ρ_{k} / 100)^{m_{k}} .$

$β_{k} = 1 + ρ_{k} / 100$ を定義し、 $r_{k}$ を債券 k のトータルリターンとして定義します。

$r_{k} = (1 + ρ_{k} / 100)^{m_{k}} = β_{k}^{m_{k}} .$

x = optimvar('x',nPtotal,'LowerBound',0);
% Total returns
r = rt.^Maturity;

目的関数

ゴールは T 年目の終わりに回収される金額を最大化する投資を選択することです。プロットから、投資は中間のさまざまな年に回収され、再投資されることがわかります。T 年目の終わりに投資 5、7、および 8 からリターンを回収でき、最終資産額は次のように表されます。

$\max_{x} x_{5} r_{5} + x_{7} r_{7} + x_{8} r_{8}$

最大化を目的とする最適化問題を作成し、目的関数を含めます。

interestprob = optimproblem('ObjectiveSense','maximize');
interestprob.Objective = x(5)*r(5) + x(7)*r(7) + x(8)*r(8);

線形制約: 保有金額を超えて投資しない

毎年、一定金額を債券の購入に利用できます。まず 1 年目には、初期資本を購入オプション $x_{1}$ および $x_{6}$ に投資できるため、次のようになります。

$x_{1} + x_{6} = {C a p i t a l}_{0}$

その後の年では、満期になった債券からリターンを回収し、利用可能な新しい債券に再投資するため、次の方程式系が得られます。

$\begin{array}{ll} x_{2} + x_{8} + x_{9} & = r_{1} x_{1} \\ x_{3} & = r_{2} x_{2} \\ x_{4} & = r_{3} x_{3} \\ x_{5} + x_{7} & = r_{4} x_{4} + r_{6} x_{6} + r_{9} x_{9} \end{array}$

investconstr = optimconstr(T,1);
investconstr(1) = x(1) + x(6) == Capital_0;
investconstr(2) = x(2) + x(8) + x(9) == r(1)*x(1);
investconstr(3) = x(3) == r(2)*x(2);
investconstr(4) = x(4) == r(3)*x(3);
investconstr(5) = x(5) + x(7) == r(4)*x(4) + r(6)*x(6) + r(9)*x(9);
interestprob.Constraints.investconstr = investconstr;

範囲制約: 借入をしない

各投資金額が正の値でなければならないため、解ベクトル $x$ の各エントリは正の値でなければなりません。解ベクトル $x$ に下限を設定して、この制約を含めます。解ベクトルには明示的な上限はありません。

x.LowerBound = 0;

問題を解く

債券に投資できる金額に制約を付けずにこの問題を解きます。このような線形計画問題を解くには、内点法アルゴリズムを使用できます。

options = optimoptions('linprog','Algorithm','interior-point');
[sol,fval,exitflag] = solve(interestprob,'options',options)

Solving problem using linprog.

Solution found during presolve.

sol = struct with fields:
    x: [9×1 double]

fval = 1.2625e+03

exitflag = 
    OptimalSolution

解の可視化

終了フラグは、ソルバーによって最適解が見つかったことを示しています。2 番目の出力引数として返される値 fval は、最終資産額に相当します。投資額の最終的な合計と時間の経過に合わせた投資配分を確認します。

fprintf('After %d years, the return for the initial $%g is $%g \n',...
    T,Capital_0,fval);

After 5 years, the return for the initial $1000 is $1262.48

plotInvestments(N,PurchaseYears,Maturity,InterestRates,sol.x)

保有債券に制限がある場合の最適な投資

投資を分散するために、1 つの債券に投資する金額を、その年の総資本 (その時点で満期になっている債券のリターンを含む) の特定の割合 Pmax に制限するよう選択することもできます。次の連立不等式が得られます。

$\begin{array}{ll} x_{1} & \leq P m a x \times {C a p i t a l}_{0} \\ x_{2} & \leq P m a x \times (β_{1} x_{1} + β_{6} x_{6}) \\ x_{3} & \leq P m a x \times (β_{2} x_{2} + β_{6}^{2} x_{6} + β_{8} x_{8} + β_{9} x_{9}) \\ x_{4} & \leq P m a x \times (β_{3} x_{3} + β_{6}^{3} x_{6} + β_{8}^{2} x_{8} + β_{9}^{2} x_{9}) \\ x_{5} & \leq P m a x \times (β_{4} x_{4} + β_{6}^{4} x_{4} + β_{8}^{3} x_{8} + β_{9}^{3} x_{9}) \\ x_{6} & \leq P m a x \times {C a p i t a l}_{0} \\ x_{7} & \leq P m a x \times (β_{4} x_{4} + β_{6}^{4} x_{4} + β_{8}^{3} x_{8} + β_{9}^{3} x_{9}) \\ x_{8} & \leq P m a x \times (β_{1} x_{1} + β_{6} x_{6}) \\ x_{9} & \leq P m a x \times (β_{1} x_{1} + β_{6} x_{6}) \end{array}$

% Maximum percentage to invest in any bond
Pmax = 0.6;

constrlimit = optimconstr(nPtotal,1);
constrlimit(1) = x(1) <= Pmax*Capital_0;
constrlimit(2) = x(2) <= Pmax*(rt(1)*x(1) + rt(6)*x(6));
constrlimit(3) = x(3) <= Pmax*(rt(2)*x(2) + rt(6)^2*x(6) + rt(8)*x(8) + rt(9)*x(9));
constrlimit(4) = x(4) <= Pmax*(rt(3)*x(3) + rt(6)^3*x(6) + rt(8)^2*x(8) + rt(9)^2*x(9));
constrlimit(5) = x(5) <= Pmax*(rt(4)*x(4) + rt(6)^4*x(6) + rt(8)^3*x(8) + rt(9)^3*x(9));
constrlimit(6) = x(6) <= Pmax*Capital_0;
constrlimit(7) = x(7) <= Pmax*(rt(4)*x(4) + rt(6)^4*x(6) + rt(8)^3*x(8) + rt(9)^3*x(9));
constrlimit(8) = x(8) <= Pmax*(rt(1)*x(1) + rt(6)*x(6));
constrlimit(9) = x(9) <= Pmax*(rt(1)*x(1) + rt(6)*x(6));

interestprob.Constraints.constrlimit = constrlimit;

どの資産への投資も 60% を超えないようにして問題を解きます。結果として得られる購入額をプロットします。この制約なしに投資する場合より最終資産額が少なくなることがわかります。

[sol,fval] = solve(interestprob,'options',options);

Solving problem using linprog.

Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in feasible directions, to within the function tolerance, and constraints are satisfied to within the constraint tolerance.

fprintf('After %d years, the return for the initial $%g is $%g \n',...
    T,Capital_0,fval);

After 5 years, the return for the initial $1000 is $1207.78

plotInvestments(N,PurchaseYears,Maturity,InterestRates,sol.x)

任意のサイズのモデル

一般的な問題のモデルを作成します。T = 30 年、利率が 1 ～ 6% でランダムに生成された 400 個の債券を使用してこのモデルを説明します。この設定の場合、決定変数が 430 個ある線形計画問題になります。

% For reproducibility
rng default 
% Initial amount of money
Capital_0 = 1000;
% Time period in years
T = 30;
% Number of bonds
N = 400;
% Total number of buying oportunities
nPtotal = N + T;
% Generate random maturity durations
Maturity = randi([1 T-1],nPtotal,1);
% Bond 1 has a maturity period of 1 year
Maturity(1:T) = 1;
% Generate random yearly interest rate for each bond
InterestRates = randi(6,nPtotal,1);
% Bond 1 has an interest rate of 0 (not invested)
InterestRates(1:T) = 0;
% Return after one year of interest
rt = 1 + InterestRates/100;
% Compute the return at the end of the maturity period for each bond:
r = rt.^Maturity;

% Generate random purchase years for each option
PurchaseYears = zeros(nPtotal,1);
% Bond 1 is available for purchase every year
PurchaseYears(1:T)=1:T;
for i=1:N
    % Generate a random year for the bond to mature before the end of
    % the T year period
    PurchaseYears(i+T) = randi([1 T-Maturity(i+T)+1]);
end

% Compute the years where each bond reaches maturity at the end of the year
MaturityYears = PurchaseYears + Maturity - 1;

債券を購入または売却できる時点を計算します。行列 buyindex には購入可能な時点が保持され、行列 sellindex には各債券の売却可能な時点が保持されます。

buyindex = false(nPtotal,T); % allocate nPtotal-by-T matrix
for ii = 1:T
    buyindex(:,ii) = PurchaseYears == ii;
end
sellindex = false(nPtotal,T);
for ii = 1:T
    sellindex(:,ii) = MaturityYears == ii;
end

債券に対応する最適化変数を設定します。

x = optimvar('x',nPtotal,1,'LowerBound',0);

最適化問題と目的関数を作成します。

interestprob = optimproblem('ObjectiveSense','maximize');
interestprob.Objective = sum(x(sellindex(:,T)).*r(sellindex(:,T)));

便宜上、一時配列 xBuy を作成します。この配列の列は、各期間で購入可能な債券を表しています。

xBuy = repmat(x,1,T).*double(buyindex);

同様に、一時配列 xSell を作成します。この配列の列は、各期間で売却可能な債券を表しています。

xSell = repmat(x,1,T).*double(sellindex);

これらの債券の売却によって生じる収益は次のとおりです。

xReturnFromSell = xSell.*repmat(r,1,T);

各期間の投資額を直前の期間に売却した金額とする制約を作成します。

interestprob.Constraints.InitialInvest = sum(xBuy(:,1)) == Capital_0;
interestprob.Constraints.InvestConstraint = sum(xBuy(:,2:T),1) == sum(xReturnFromSell(:,1:T-1),1);

保有債券に制限がない場合の解

問題を解きます。

tic
[sol,fval,exitflag] = solve(interestprob,'options',options);

Solving problem using linprog.

Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in feasible directions, to within the function tolerance, and constraints are satisfied to within the constraint tolerance.

toc

Elapsed time is 0.033693 seconds.

投資の成果はどうでしょうか。

fprintf('After %d years, the return for the initial $%g is $%g \n',...
    T,Capital_0,fval);

After 30 years, the return for the initial $1000 is $5167.58

保有債券に制限がある場合の解

各資産への投資の割合を制限する制約を作成するには、各時点でアクティブな債券を追跡する行列を設定します。各投資額が合計値の Pmax 倍未満でなければならないという制約を表すには、各時点でのそれぞれの投資の値を追跡する行列を設定します。このより大きな問題では、保持できる最大の割合を 0.4 に設定します。

Pmax = 0.4;

債券を保持できる時点に対応する行列 active およびアクティブな債券ごとの累積期間を保持する行列 cactive を作成します。時間 t での債券 j の値は x(j)*(rt^cactive) になります。

active = double(buyindex | sellindex);
for ii = 1:T
    active(:,ii) = double((ii >= PurchaseYears) & (ii <= MaturityYears));
end
cactive = cumsum(active,2);
cactive = cactive.*active;

エントリ (j,p) が期間 p の債券 j の値を表す行列を作成します。

bondValue = repmat(x, 1, T).*active.*(rt.^(cactive));

制限がある保有債券に制約を課すことができるように各期間の投資の合計値を求めます。mvalue は各期間の終わりにおいてすべての債券に投資されている金額で、nPtotal 行 T 列の行列です。moneyavailable は、その期間の始まりにおいて投資されている金額の債券全体に対する合計、つまり、各時点のポートフォリオの価値です。

constrlimit = optimconstr(nPtotal,T);
constrlimit(:,1) = xBuy(:,1) <= Pmax*Capital_0;
constrlimit(:,2:T) = xBuy(:,2:T)  <= repmat(Pmax*sum(bondValue(:,1:T-1),1), nPtotal, 1).*double(buyindex(:,2:T));
interestprob.Constraints.constrlimit = constrlimit;

保有債券に制限がある問題を解きます。

tic
[sol,fval,exitflag] = solve(interestprob,'options',options);

Solving problem using linprog.

Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in feasible directions, to within the function tolerance, and constraints are satisfied to within the constraint tolerance.

toc

Elapsed time is 0.627319 seconds.

fprintf('After %d years, the return for the initial $%g is $%g \n',...
    T,Capital_0,fval);

After 30 years, the return for the initial $1000 is $5095.26

ソルバーのスピードを向上させるために、双対シンプレックス法アルゴリズムを使用してみます。

options = optimoptions('linprog','Algorithm','dual-simplex');
tic
[sol,fval,exitflag] = solve(interestprob,'options',options);

Solving problem using linprog.

Optimal solution found.

toc

Elapsed time is 0.149014 seconds.

fprintf('After %d years, the return for the initial $%g is $%g \n',...
    T,Capital_0,fval);

After 30 years, the return for the initial $1000 is $5095.26

このケースでは、双対シンプレックス法アルゴリズムは同じ解を得るまでの時間を短縮しました。

定性結果分析

解の特性を確認するため、この解を、初期資金の全額を利率 6% (最大利率) の 1 つの債券に 30 年間にわたって投資した場合に得られる金額 fmax と比較します。最終資産額に対応する同等の利率を計算することもできます。

% Maximum amount
fmax = Capital_0*(1+6/100)^T;
% Ratio (in percent)
rat = fval/fmax*100;
% Equivalent interest rate (in percent)
rsol = ((fval/Capital_0)^(1/T)-1)*100;

fprintf(['The amount collected is %g%% of the maximum amount $%g '...
    'that you would obtain from investing in one bond.\n'...
    'Your final wealth corresponds to a %g%% interest rate over the %d year '...
    'period.\n'], rat, fmax, rsol, T)

The amount collected is 88.7137% of the maximum amount $5743.49 that you would obtain from investing in one bond.
Your final wealth corresponds to a 5.57771% interest rate over the 30 year period.

plotInvestments(N,PurchaseYears,Maturity,InterestRates,sol.x,false)