Максимизируйте долгосрочные инвестиции, используя линейное программирование: основанное на решателе

В этом примере показано, как использовать linprog решатель в Optimization Toolbox ®, чтобы решить инвестиционную задачу с детерминированными возвратами в течение фиксированного количества лет T. Задача состоит в том, чтобы распределить ваши деньги над доступными инвестициями, чтобы максимизировать ваше окончательное богатство. Этот пример использует подход, основанный на решателе.

Формулировка задачи

Предположим, что у вас есть начальная сумма денег Capital_0 инвестировать в течение определенного периода времени T лет в N облигации с нулевым купоном. Каждая облигация выплачивает процентную ставку, которая соединяется каждый год, и выплачивает основной плюс совокупный процент в конце периода погашения. Цель состоит в том, чтобы максимизировать общую сумму денег после T лет.

Вы можете включить ограничение, что ни одна инвестиция не больше определенной части вашего общего капитала.

В этом примере сначала показана настройка задачи в небольшом случае, а затем формулируется общий случай.

Можно смоделировать это как задачу линейного программирования. Поэтому, чтобы оптимизировать свое богатство, сформулируйте задачу для решения linprog решатель.

Вводный пример

Начните с небольшого примера:

  • Стартовая сумма для инвестирования Capital_0 $1000.

  • Период времени T составляет 5 лет.

  • Количество облигаций N равен 4.

  • Чтобы смоделировать неинвестированные деньги, имейте одну опцию B0 доступный каждый год, который имеет период погашения 1 год и процентную ставку 0%.

  • Облигация 1, обозначенная B1, может быть приобретена в 1 году, имеет срок погашения 4 года и процентную ставку 2%.

  • Облигация 2, обозначенная B2, может быть приобретена в 5 году, имеет срок погашения 1 год и процентную ставку 4%.

  • Облигация 3, обозначенная B3, может быть приобретена во 2 году, имеет срок погашения 4 года и процентную ставку 6%.

  • Облигация 4, обозначенная B4, может быть приобретена в году 2, имеет срок погашения 3 года и процентную ставку 6%.

Путем разделения первого B0 опции на 5 облигаций с периодом погашения 1 год и процентной ставкой 0%, эта проблема может быть эквивалентно смоделирована как имеющая в общей сложности 9 доступных облигаций, таких что для k=1..9

  • Введите k от вектора PurchaseYears представляет год, в который k облигаций доступно для покупки.

  • Введите k от вектора Maturity представляет период погашения mk облигационных k.

  • Введите k от вектора InterestRates представляет процентную ставку ρ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];
% Interest rates
InterestRates = [0;0;0;0;0;2;4;6;6];

plotInvestments(N,PurchaseYears,Maturity,InterestRates)

Переменные принятия решений

Представьте свои переменные принятия решений вектором x, где x(k) - долларовая сумма инвестиций в облигации k, для k=1..9. После погашения, выплата за инвестицию x(k) является

x(k)(1+ρk/100)mk.

Определить rk как общий возврат облигаций k:

rk=(1+ρk/100)mk.

% Total returns
finalReturns = (1+InterestRates/100).^Maturity;

Целевая функция

Цель состоит в том, чтобы выбрать инвестиции, чтобы максимизировать количество денег, собранных в конце года T. Из графика видно, что инвестиции собираются в различные промежуточные годы и реинвестируются. В конце года T, деньги, возвращенные от инвестиций 5, 7 и 8, могут быть собраны и представляют ваше окончательное богатство:

maxxx5r5+x7r7+x8r8

Чтобы поместить эту задачу в форму linprog решает, превращает эту задачу максимизации в задачу минимизации с помощью негатива коэффициентов x(j):

minxfTx

с

f=[0;0;0;0;-r5;0;-r7;-r8;0]

f = zeros(nPtotal,1);
f([5,7,8]) = [-finalReturns(5),-finalReturns(7),-finalReturns(8)];

Линейные ограничения: Инвестируйте не больше, чем у вас

Каждый год у вас есть определенная сумма денег для приобретения облигаций. Начиная с 1 года, вы можете инвестировать начальный капитал в опции покупки x1 и x6, так:

x1+x6=Capital0

Затем в течение следующих лет вы собираете возвраты от погашения облигаций и реинвестируете их в новые доступные облигации, чтобы получить систему уравнений:

x2+x8+x9=r1x1x3=r2x2x4=r3x3x5+x7=r4x4+r6x6+r9x9

Запишите эти уравнения в форму Aeqx=beq, где каждая строка Aeq матрица соответствует равенству, которое должно быть удовлетворено в том году:

Aeq=[100001000-r1100000110-r2100000000-r3100000000-r41-r610-r9]

beq=[Capital0000]

Aeq = spalloc(N+1,nPtotal,15);
Aeq(1,[1,6]) = 1;
Aeq(2,[1,2,8,9]) = [-1,1,1,1];
Aeq(3,[2,3]) = [-1,1];
Aeq(4,[3,4]) = [-1,1];
Aeq(5,[4:7,9]) = [-finalReturns(4),1,-finalReturns(6),1,-finalReturns(9)];

beq = zeros(T,1);
beq(1) = Capital_0;

Связанные ограничения: нет заимствований

Поскольку каждая инвестированная сумма должна быть положительной, каждая запись в векторе решения x должен быть положительным. Включите это ограничение путем установки нижней границы lb на векторе решения x. Явная верхняя граница вектора решения отсутствует. Таким образом установите верхнюю границу ub в пустой.

lb = zeros(size(f));
ub = [];

Решите задачу

Решите эту проблему без ограничений на сумму, которую вы можете инвестировать в облигацию. Алгоритм внутренней точки может использоваться, чтобы решить этот тип задачи линейного программирования.

options = optimoptions('linprog','Algorithm','interior-point');
[xsol,fval,exitflag] = linprog(f,[],[],Aeq,beq,lb,ub,options);
Solution found during presolve.

Визуализация решения

Выходной флаг равен 1, что указывает на то, что решатель нашел решение. Значение -fval, вернулся как второй linprog выходной аргумент, соответствует конечному богатству. Стройте графики инвестиций с течением времени.

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,xsol)

Оптимальные инвестиции с ограниченными холдингами

Чтобы диверсифицировать свои инвестиции, вы можете принять решение ограничить сумму вложений в любую одну облигацию определенным процентом Pmax от общего капитала в том году (включая возвраты по облигациям, которые в настоящее время находятся в периоде погашения). Вы получаете следующую систему неравенств:

x1Pmax×Capital0x2Pmax×(ρ1x1+ρ6x6)x3Pmax×(ρ2x2+ρ62x6+ρ8x8+ρ9x9)x4Pmax×(ρ3x3+ρ63x6+ρ82x8+ρ92x9)x5Pmax×(ρ4x4+ρ64x4+ρ83x8+ρ93x9)x6Pmax×Capital0x7Pmax×(ρ4x4+ρ64x4+ρ83x8+ρ93x9)x8Pmax×(ρ1x1+ρ6x6)x9Pmax×(ρ1x1+ρ6x6)

Поместите эти неравенства в матричный вид Ax <= b.

Чтобы настроить систему неравенств, сначала сгенерируйте матрицу yearlyReturns который содержит возврат для облигации, индексируемой i в году j в строке i и столбце j. Представьте эту систему как разреженную матрицу.

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

% Build the return for each bond over the maturity period as a sparse
% matrix
cumMaturity = [0;cumsum(Maturity)];
xr = zeros(cumMaturity(end-1),1);
yr = zeros(cumMaturity(end-1),1);
cr = zeros(cumMaturity(end-1),1);
for i = 1:nPtotal
    mi = Maturity(i); % maturity of bond i
    pi = PurchaseYears(i); % purchase year of bond i
    idx = cumMaturity(i)+1:cumMaturity(i+1); % index into xr, yr and cr
    xr(idx) = i; % bond index
    yr(idx) = pi+1:pi+mi; % maturing years
    cr(idx) = (1+InterestRates(i)/100).^(1:mi); % returns over the maturity period
end
yearlyReturns = sparse(xr,yr,cr,nPtotal,T+1);

% Build the system of inequality constraints
A = -Pmax*yearlyReturns(:,PurchaseYears)'+ speye(nPtotal);

% Left-hand side
b = zeros(nPtotal,1);
b(PurchaseYears == 1) = Pmax*Capital_0;

Решить проблему путем инвестирования не более 60% в любой один актив. Постройте график полученных покупок. Заметьте, что ваше конечное богатство меньше, чем инвестиции без этого ограничения.

[xsol,fval,exitflag] = linprog(f,A,b,Aeq,beq,lb,ub,options);
Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in
feasible directions, to within the selected value of the function tolerance,
and constraints are satisfied to within the selected value of 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,xsol)

Модель произвольного размера

Создайте модель для общей версии задачи. Проиллюстрировать его используя T = 30 лет и 400 случайным образом сгенерированных облигаций с процентными ставками от 1 до 6%. Эта настройка приводит к задаче линейного программирования с 430 переменными принятия решений. Система ограничений равенства представлена разреженной матрицей Aeq размерности 30 на 430, и система неравенств представлена разреженной матрицей A Размерность 430 на 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;
% Compute the return at the end of the maturity period for each bond:
finalReturns = (1+InterestRates/100).^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
SaleYears = PurchaseYears + Maturity;

% Initialize f to 0
f = zeros(nPtotal,1);
% Indices of the sale oportunities at the end of year T
SalesTidx = SaleYears==T+1;
% Expected return for the sale oportunities at the end of year T
ReturnsT = finalReturns(SalesTidx);
% Objective function
f(SalesTidx) = -ReturnsT;


% Generate the system of equality constraints.
% For each purchase option, put a coefficient of 1 in the row corresponding
% to the year for the purchase option and the column corresponding to the
% index of the purchase oportunity
xeq1 = PurchaseYears;
yeq1 = (1:nPtotal)';
ceq1 = ones(nPtotal,1);

% For each sale option, put -\rho_k, where \rho_k is the interest rate for the
% associated bond that is being sold, in the row corresponding to the
% year for the sale option and the column corresponding to the purchase
% oportunity
xeq2 = SaleYears(~SalesTidx);
yeq2 = find(~SalesTidx);
ceq2 = -finalReturns(~SalesTidx);

% Generate the sparse equality matrix
Aeq = sparse([xeq1; xeq2], [yeq1; yeq2], [ceq1; ceq2], T, nPtotal);

% Generate the right-hand side
beq = zeros(T,1);
beq(1) = Capital_0;

% Build the system of inequality constraints
% Maximum percentage to invest in any bond
Pmax = 0.4;

% Build the returns for each bond over the maturity period
cumMaturity = [0;cumsum(Maturity)];
xr = zeros(cumMaturity(end-1),1);
yr = zeros(cumMaturity(end-1),1);
cr = zeros(cumMaturity(end-1),1);
for i = 1:nPtotal
    mi = Maturity(i); % maturity of bond i
    pi = PurchaseYears(i); % purchase year of bond i
    idx = cumMaturity(i)+1:cumMaturity(i+1); % index into xr, yr and cr
    xr(idx) = i; % bond index
    yr(idx) = pi+1:pi+mi; % maturing years
    cr(idx) = (1+InterestRates(i)/100).^(1:mi); % returns over the maturity period
end
yearlyReturns = sparse(xr,yr,cr,nPtotal,T+1);

% Build the system of inequality constraints
A = -Pmax*yearlyReturns(:,PurchaseYears)'+ speye(nPtotal);

% Left-hand side
b = zeros(nPtotal,1);
b(PurchaseYears==1) = Pmax*Capital_0;

% Add the lower-bound constraints to the problem.
lb = zeros(nPtotal,1);

Решение без предела удержания

Во-первых, решите задачу линейного программирования без ограничений неравенства, используя алгоритм внутренней точки.

% Solve the problem without inequality constraints
options = optimoptions('linprog','Algorithm','interior-point');
tic
[xsol,fval,exitflag] = linprog(f,[],[],Aeq,beq,lb,[],options);
Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in
feasible directions, to within the selected value of the function tolerance,
and constraints are satisfied to within the selected value of the constraint
tolerance.
toc
Elapsed time is 0.053286 seconds.
fprintf('\nAfter %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 

Решение с ограниченными холдингами

Теперь решите проблему с ограничениями неравенства.

% Solve the problem with inequality constraints
options = optimoptions('linprog','Algorithm','interior-point');
tic
[xsol,fval,exitflag] = linprog(f,A,b,Aeq,beq,lb,[],options);
Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in
feasible directions, to within the selected value of the function tolerance,
and constraints are satisfied to within the selected value of the constraint
tolerance.
toc
Elapsed time is 1.360809 seconds.
fprintf('\nAfter %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 

Даже при том, что количество ограничений увеличилось на порядок 10, время, которое решатель должен найти решение, увеличилось на порядок 100. Это расхождение эффективности частично вызвано плотными столбцами в системе неравенства, показанной в матрице A. Эти столбцы соответствуют связям с длительным периодом погашения, как показано на следующем графике.

% Number of nonzero elements per column
nnzCol = sum(spones(A));

% Plot the maturity length vs. the number of nonzero elements for each bond
figure;
plot(Maturity,nnzCol,'o');
xlabel('Maturity period of bond k')
ylabel('Number of nonzero in column k of A')

Figure contains an axes. The axes contains an object of type line.

Плотные столбцы в ограничениях приводят к плотным блокам во внутренних матрицах решателя, что приводит к потере эффективности его разреженных методов. Чтобы ускорить решатель, попробуйте алгоритм двойного симплекса, который менее чувствителен к плотности столбца.

% Solve the problem with inequality constraints using dual simplex
options = optimoptions('linprog','Algorithm','dual-simplex');
tic
[xsol,fval,exitflag] = linprog(f,A,b,Aeq,beq,lb,[],options);
Optimal solution found.
toc
Elapsed time is 0.244413 seconds.
fprintf('\nAfter %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 

В этом случае алгоритму dual-simplex потребовалось гораздо меньше времени, чтобы получить одно и то же решение.

Качественный анализ результатов

Чтобы почувствовать решение, найденное linprog, сравните его с суммой fmax что вы получите, если сможете инвестировать все свои стартовые деньги в одну облигацию с процентной ставкой 6% (максимальная процентная ставка) за полный 30-летний период. Можно также вычислить эквивалентную процентную ставку, соответствующую вашему конечному богатству.

% 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,xsol,false)

Похожие темы