Рассчитать доходность облигаций к погашению (YTM) в C ++

Question

Я хочу рассчитать доходность облигации к погашению с учетом цены, используя либо бисекцию, либо секущий метод. Я знаю, что для этого есть рецепты на С ++, но я не могу понять, что не так с моим собственным кодом. Оба метода создают бесконечное число l oop. Когда я пытался перебирать переменные, итерационные решения доходности к погашению увеличивались до 80% и выше.

#include <iostream>
#include <math.h>

using namespace std;

class bond {
private:
    double principal;
    double coupon;
    double timeToMaturity;

public:
    bond(double principal, double coupon, double timeToMaturity) {
        this->principal = principal;
        this->coupon = coupon;
        this->timeToMaturity = timeToMaturity;
    }

    double getprice(double YTM);
    double getytm(double price);
    double getytmsec(double price);
};

bool isZero(double x)
{
    if (fabs(x) < 1e-10)
        return true;
    return false;
}

double bond::getprice(double YTM) {
    double pvPrincipal;
    double pvCoupon;
    double factor = 0;

    pvPrincipal = this->principal / pow(1 + YTM, this->timeToMaturity);

    for (int i = 0; (this->timeToMaturity - i) > 0; i++) {
        double denom = pow(1 + YTM, this->timeToMaturity - i);
        factor += 1 / denom;
    }

    pvCoupon = this->coupon * factor;

    return pvPrincipal + pvCoupon;
}

// Bisection method
double bond::getytm(double price) {
    double low = 0;
    double high = 1;

    double f0 = getprice(low);
    double f2 = 1;
    double x2 = 0;

    while (!isZero(f2)) {
        x2 = (low + high) / 2;
        f2 = getprice(x2) - price;
        if (f2 < 0) {
            low = x2;
        }
        else {
            high = x2;
        }
    }

    return x2;
}

// Secant method
double bond::getytmsec(double price) {
    double x1 = price;
    double x2 = price + 0.25;
    double f1 = this->getprice(x1);
    double f2 = this->getprice(x2);
    for (; !isZero(f1 - price); x2 = x1, x1 = price) {
        f1 = getprice(x1);
        f2 = getprice(x2);
        price = x1 - f1 * (x1 - x2) / (f1 - f2);
    }
    return price;
}

int main() {
    bond myBond = { 1000, 25, 6 };
    cout << "YTM is " << myBond.getytm(950) << endl;
    cout << "YTM is " << myBond.getytmsec(950) << endl;

    return 0;
}

Answer 1

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

Проблема в том, что нужно найти ноль для функции f(x) = getprice(x) - price.

В целом метод деления пополам: начинать с интервала [low, high], где f(low) и f(high) имеют разные знаки (один не положительный, другой неотрицательный). Это означает, что он содержит ноль. Затем выберите левый или правый подинтервал на основе значения функции в средней точке, чтобы сохранить это свойство.

В этом случае функция однотонная c и не увеличивается, поэтому мы знаем, что f(low) должно быть большим (неотрицательным) числом, а f(high) должно быть меньшим (неотрицательным) числом. Поэтому мы должны выбрать левый подинтервал, если f(midpoint) отрицательно, и выбрать правый подинтервал, если f(midpoint) положительно. Но код делает противоположное, выбирая правильный подинтервал, если f(midpoint) отрицателен:

    x2 = (low + high) / 2;
    f2 = getprice(x2) - price;
    if (f2 < 0) {
        low = x2;
    }

Таким образом, вы выбираете все меньшие и меньшие правые подинтервалы с в конечном итоге [low, high] = [1, 1], и это бесконечный l oop. Замените f2 < 0 на f2 > 0.

. В секущем методе обычно используются две нулевые "оценки" x_k и x_{k-1} и используется повторение для поиска лучшей "оценки" x_{k+1}. Повторение по существу использует линию между (x_{k-1}, f(x_{k-1}) и (x_k, f(x_k)) и ищет, где эта линия пересекает ноль.

В предоставленном коде есть несколько проблем. Во-первых, на важном этапе:

    price = x1 - f1 * (x1 - x2) / (f1 - f2);

, где x1 и x2 - текущие и предыдущие оценки, f1 - getprice(x1) и f2 - getprice(x2). Важно отметить, что f1 - это не f(x1), где f - это функция, ноль которой мы хотим. Это не секущая формула. Первая часть второго члена должна быть значением функции в x1, то есть f1 - price, а не f1:

    ... = x1 - (f1 - price) * (x1 - x2) / (f1 - f2);

Во-вторых, вы присваиваете это price и тем самым теряя фактическое значение price, которое вам нужно на каждой итерации.

В-третьих, начальные предположения о доходности: price и price + 0.25. Они настолько далеки от фактического значения, что это становится проблемой (ноль - это доходность между 0 и 1). Попробуйте 0 и 1.

Многое из этого можно избежать, не смешивая многие проблемы. Вы можете вычленить логику c нахождения нуля функции из фактического идентификатора функции. Например, один шаг деления пополам:

template<typename Function>
constexpr auto bisection_step(double low, double high, Function f) -> std::pair<double, double>
{
    assert(std::isfinite(high));
    assert(std::isfinite(low));
    assert(low < high);
    assert(f(low) * f(high) <= 0.);

    auto mid = midpoint(low, high);
    if (f(low) * f(mid) <= 0.)
        return {low, mid};
    else
        return {mid, high};
}

Это позволяет вам указать предположения в форме утверждений или проверок, которые генерируют исключения или возвращают коды ошибок. Это также делает логику c более понятной, поэтому менее вероятно, что вы выберете неправильный подинтервал. И даже если бы кто-то сделал, утверждение сработало бы.