Нелинейная логистическая регрессия

Открыть скрипт

Этот пример показывает два способа подбирать нелинейную модель логистической регрессии. Первый метод использует наибольшее правдоподобие (ML), и второй метод использует обобщенные наименьшие квадраты (GLS) через функциональный fitnlm от Statistics and Machine Learning Toolbox(TM).

Описание проблемы

Логистическая регрессия является специальным типом регрессии, в которой цель состоит в том, чтобы смоделировать вероятность чего-то как функция других переменных. Рассмотрите набор векторов предиктора $x_1,\ldots,x_N$ , где $N$ количество наблюдений и $x_i$ вектор-столбца, содержит значения $d$ предикторов для $i$ th наблюдения. Переменная отклика для $x_i$ - $Z_i$ то, где $Z_i$ представляет Биномиальную случайную переменную параметрами $n$ , количеством испытаний, и $\mu_i$ , вероятность успеха для испытания $i$ . Нормированная переменная отклика $Y_i = Z_i/n$ - пропорция успехов в $n$ испытаниях за наблюдение $i$ . Примите, что ответы $Y_i$ независимы для $i = 1,\ldots,N$ . Для каждого $i$ :

$E(Y_i) = \mu_i$

$Var(Y_i) = \frac{\mu_i(1-\mu_i)}{n}$ .

Рассмотрите моделирование $\mu_i$ как функцию переменных предикторов $x_i$ .

В линейной логистической регрессии можно использовать функциональный fitglm к модели $\mu_i$ как функция $x_i$ можно следующим образом:

$\log{\left( \frac{\mu_i}{1 - \mu_i} \right)} = x_i^T \beta$

с $\beta$ представлением набора коэффициентов, умножающих предикторы в $x_i$ . Однако предположите, что вам нужна нелинейная функция на правой стороне:

$\log{\left( \frac{\mu_i}{1 - \mu_i} \right)} = f(x_i,\beta).$

Существуют функции в Statistics and Machine Learning Toolbox(TM) для того, чтобы подбирать нелинейные модели регрессии, но не для того, чтобы подбирать нелинейные модели логистической регрессии. В этом примере показано, как можно использовать функции тулбокса, чтобы подбирать те модели.

Прямое наибольшее правдоподобие (ML)

Подход ML максимизирует логарифмическую вероятность наблюдаемых данных. Вероятность легко вычисляется с помощью Биномиальной вероятности (или плотность) функция, как вычислено binopdf функция.

Обобщенные наименьшие квадраты (GLS)

Можно оценить нелинейную модель логистической регрессии использование функционального fitnlm. Эта сила кажется удивительной сначала начиная с fitnlm не вмещает Биномиальное распределение или любые функции ссылки. Однако fitnlm может использовать Обобщенные наименьшие квадраты (GLS) для оценки модели, если вы задаете среднее значение и отклонение ответа. Если GLS сходится, то он решает тот же набор нелинейных уравнений для оценки $\beta$ , как решено ML. Можно также использовать GLS для оценки квазивероятности обобщенных линейных моделей. Другими словами, мы должны получить те же или эквивалентные решения от GLS и ML. Чтобы реализовать оценку GLS, обеспечьте нелинейную функцию, чтобы соответствовать, и функция отклонения для Биномиального распределения.

Средняя или функция модели

Функция модели описывает как $\mu_i$ изменения с $\beta$ . Для fitnlm, функция модели:

$\mu_i = \frac{1}{1 \,\, + \,\, \mbox{exp}\left\{-f(x_i,\beta)\right\}}$

Функция веса

fitnlm принимает веса наблюдения как указатель на функцию с помощью 'Weights' аргумент пары "имя-значение". При использовании этой опции, fitnlm принимает следующую модель:

$E(Y_i) = \mu_i$

$Var(Y_i) = \frac{\sigma^2}{w(\mu_i)}$

где ответы $Y_i$ приняты, чтобы быть независимыми, и $w$ пользовательский указатель на функцию, который принимает $\mu_i$ и возвращает вес наблюдения. Другими словами, веса обратно пропорциональны отклонению ответа. Для Биномиального распределения, используемого в модели логистической регрессии, создайте функцию веса можно следующим образом:

$w(\mu_i) = \frac{1}{Var(y_i)} = \frac{n}{\mu_i (1 - \mu_i)}$

fitnlm моделирует отклонение ответа $Y_i$ как $\sigma^2/w(\mu_i)$ , где $\sigma^2$ дополнительный параметр, который присутствует по оценке GLS, но отсутствующий в модели логистической регрессии. Однако это обычно не влияет на оценку $\beta$ , и она обеспечивает "дисперсионный параметр", чтобы проверять при условии, что $Z_i$ значения имеют Биномиальное распределение.

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

Сгенерируйте данные в качестве примера

Чтобы проиллюстрировать различия между ML и подбором кривой GLS, сгенерируйте некоторые данные в качестве примера. Примите, что это $x_i$ является одномерным, и предположите, что истинная функция $f$ в нелинейной модели логистической регрессии является моделью Michaelis-Menten, параметрированной $2 \times 1$ вектором $\beta$ :

$f(x_i,\beta) = \frac{\beta_1 x_i}{\beta_2 + x_i}.$

    myf = @(beta,x) beta(1)*x./(beta(2) + x);

Создайте функцию модели, которая задает отношение между $\mu_i$ и $\beta$ .

    mymodelfun = @(beta,x) 1./(1 + exp(-myf(beta,x)));

Создайте вектор из одномерных предикторов и истинного вектора коэффициентов $\beta$ .

    rng(300,'twister');
    x    = linspace(-1,1,200)';
    beta = [10;2];

Вычислите вектор из $\mu_i$ значений для каждого предиктора.

    mu = mymodelfun(beta,x);

Сгенерируйте ответы $z_i$ от Биномиального распределения с вероятностями успеха $\mu_i$ и количеством испытаний $n$ .

    n = 50;
    z = binornd(n,mu);

Нормируйте ответы.

    y = z./n;

Подход ML

Подход ML задает отрицательную логарифмическую вероятность как функцию $\beta$ вектора, и затем минимизирует его с оптимизационной функцией, такой как fminsearch. Задайте beta0 как начальное значение для $\beta$ .

    mynegloglik = @(beta) -sum(log(binopdf(z,n,mymodelfun(beta,x))));
    beta0 = [3;3];
    opts = optimset('fminsearch');
    opts.MaxFunEvals = Inf;
    opts.MaxIter = 10000;
    betaHatML = fminsearch(mynegloglik,beta0,opts)

betaHatML =

    9.9259
    1.9720

Предполагаемые коэффициенты в betaHatML близко к истинным значениям [10;2].

Подход GLS

Подход GLS создает функцию веса для fitnlm ранее описанный.

    wfun = @(xx) n./(xx.*(1-xx));

Вызовите fitnlm с пользовательским средним значением и функциями веса. Задайте beta0 как начальное значение для $\beta$ .

    nlm = fitnlm(x,y,mymodelfun,beta0,'Weights',wfun)

nlm = 


Nonlinear regression model:
    y ~ F(beta,x)

Estimated Coefficients:
          Estimate      SE       tStat       pValue  
          ________    _______    ______    __________

    b1     9.926      0.83135     11.94     4.193e-25
    b2     1.972      0.16438    11.996    2.8182e-25


Number of observations: 200, Error degrees of freedom: 198
Root Mean Squared Error: 1.16
R-Squared: 0.995,  Adjusted R-Squared 0.995
F-statistic vs. zero model: 1.88e+04, p-value = 2.04e-226

Получите оценку $\beta$ от подходящего NonLinearModel объект nlm.

    betaHatGLS = nlm.Coefficients.Estimate

betaHatGLS =

    9.9260
    1.9720

Как в методе ML, предполагаемых коэффициентах в betaHatGLS близко к истинным значениям [10;2]. Маленькие p-значения для обоих $\beta_1$ и $\beta_2$ указывают, что оба коэффициента существенно отличаются от $0$ .

Сравните подходы GLS и ML

Сравните оценки $\beta$ .

    max(abs(betaHatML - betaHatGLS))

ans =

   1.1461e-05

Сравните адаптированные значения с помощью ML и GLS

    yHatML  = mymodelfun(betaHatML ,x);
    yHatGLS = mymodelfun(betaHatGLS,x);
    max(abs(yHatML - yHatGLS))

ans =

   1.2747e-07

ML и подходы GLS производят аналогичные решения.

Постройте адаптированные значения с помощью ML и GLS

    h = figure;
    plot(x,y,'g','LineWidth',1);
    hold on;
    plot(x,yHatML ,'b'  ,'LineWidth',1);
    plot(x,yHatGLS,'m--','LineWidth',1);
    legend('Data','ML','GLS','Location','Best');
    xlabel('x');
    ylabel('y and fitted values');
    title('Data y along with ML and GLS fits.');
    snapnow;
    close(h);

ML и GLS производят подобные подходящие значения.

Постройте оцененную нелинейную функцию с помощью ML и GLS

Постройте истинную модель для $f(x_i,\beta)$ . Добавьте график для первоначальной оценки $f(x_i,\beta)$ использования $\beta = \beta_0$ и графики для ML, и GLS основывал оценки $f(x_i,\beta)$ .

    h = figure;
    plot(x,myf(beta,x),'r','LineWidth',1);
    hold on;
    plot(x,myf(beta0,x)     ,'k'  ,'LineWidth',1);
    plot(x,myf(betaHatML ,x),'c--','LineWidth',1);
    plot(x,myf(betaHatGLS,x),'b-.','LineWidth',1);
    legend('True f','Initial f','Estimated f with ML','Estimated f with GLS','Location','Best');
    xlabel('x');
    ylabel('True and estimated f');
    title('Comparison of true f with estimated f using ML and GLS.');
    snapnow;
    close(h);

Предполагаемая нелинейная функция $f$ с помощью и ML и методов GLS близко к истинной нелинейной функции $f$ . Можно использовать подобный метод, чтобы подбирать другие нелинейные обобщенные линейные модели как нелинейная регрессия Пуассона.

Документация