векторизованная функция fminsearch - PullRequest
0 голосов
/ 26 февраля 2012

Большое спасибо, сначала. Первый код работает хорошо:
функция f = малибу (м1, м2, м3, к, л, т) F = (m1 * к + (m1 + m2) * л..) * ехр (-m3 * т.). конец

function loglik= modelmalibu(p)
global k l t x n m2 m3;
f =malibu(p,m2,m3,k,l,t);
if f==0;
    loglik0=0;
else
loglik0=(x.*log(f)+(n-x).*log(1-f));%minus likelihood
end
loglik=sum(-loglik0);
end

clear all;
global n t x k l m2 m3;
m1=0.1;
m2=0.2;
m3=0.3;
t=[1 3 6 9 12 18]';
k=[1 1 2 3 3 4]';
l=[0 0 1 3 4 5]';
y=meltem(m1,m2,m3,k,l,t);
n=100;%trial
x=y.*n;%correct replies
pstart=0.3;
[p1,modelvalue]=fminsearch(@modelmalibu,pstart);

Но подобный код для большей переменной выдает ошибку.

function w=anemon(m1,m2,m3,X,Y,k,l)
w=(m1.*k+(m1+m2).*l)+X.*exp(-m3.*Y);
end

function loglik= modelanemon(p)
global n x m2 m3 X Y k l ;
f =anemon(p,m2,m3,X,Y,k,l);
if f==0;
    loglik0=0;
else
loglik0=(x*log(f)+(n-x)*log(1-f));%minus likelihood
end
loglik=sum(-loglik0);
end

clear;
global n x Ydata kdata ldata m1 m2 m3;
%parameters
m1=0.002;
m2=0.0001;
m3=7;
%given data
Xdata=[1 3 6 9 10 12]';
Ydata=[11 13 41 81 121 181]';
kdata=[1 1 2 4 5 4]';
ldata=[1 1 3 3 4 5]';
y=anemon(m1,m2,m3,Xdata,Ydata,kdata,ldata);
n=10;
x=y.*n;
pstart=2;
[pbest,modelvalue]=fminsearch(@modelanemon,pstart);

Я на самом деле пытался использовать ваши советы, но если бы я написал неравенство вместо f == 0, первый код также упал бы.

1 Ответ

1 голос
/ 26 февраля 2012

Я думаю, что вы смешиваете две концепции.

векторизация: описывает, как написать код, чтобы он мог использовать некоторые функции ускорения в вашем процессоре.это не имеет ничего общего с fminsearch.См .: http://en.wikipedia.org/wiki/Vectorization_(parallel_computing)

Полагаю, вы хотите написать свою функцию так, чтобы она принимала вектор в качестве входных данных.Самый простой способ - просто использовать дескриптор функции следующим образом:

fh = @(x) my_complicated_function(const1, const2, x(1), x(2), x(3) )

В этом случае my_complicated_function имеет 5 входов, и вы берете первые 2 константы и вводите 3 dim вектора для другого 3. Fminsearch будет работатьс этим.

Вы бы позвонили

 x_opt = fminsearch(fh, [1,2,3])

Помимо некоторых советов по коду:

  • Не используйте == для сравнения чисел - идинапримерabs(x1-x2)<0.1 вместо
  • interpanemon выглядит очень странно - он не выполняет то, что называется интерполяцией, и, если вы посмотрите - в каждой итерации p пересчитывается, поэтому вступают в силу только последние вычисления.Поскольку это выглядит, это должно вывести постоянное значение - бесполезно оптимизировать это.
  • Использование предварительно рассчитанного Z, вероятно, не то, что вы хотите.Это уже определяет ваш оптимум.Если вы используете линейную оптимизацию, оптимум должен быть уже в Z - нет необходимости выполнять интерполяцию.
  • Вызов в вашем случае может выглядеть следующим образом:

    min_p = fminsearch (@ (x)interpanemon (5, Ydata, x, m2, m3, Z, X, Y, kdata, ldata), 1)

В целом это выглядит очень необычно - это может действительно помочь, если выобъясните идеи ПОЧЕМУ вы решили сделать это так.Также может помочь, если вы переформулируете проблему в более простой форме.

...