How to find maximum from fminunc?

5 Ansichten (letzte 30 Tage)
Boonanan Lechayakittikorn
Boonanan Lechayakittikorn am 30 Jan. 2018
Bearbeitet: Walter Roberson am 30 Jan. 2018
Hello everyone, i have problem about to find parameter of beta(0) and beta(1). i have data from generate. My problem is when i input value for generate X and y (value input are stdevu stdevex etc.) when i run code by stdevu = 0.1 and stdevex = 0.1 etc. output likely stdevu = 0.1 and stdevex = 0.3 but i change stdevu = 0.1 and stdevex = 0.1 etc. output unlikely stdevu = 0.3 and stdevex = 0.1 etc. please help me. My code is
function [EX,X,y,W,z,lam,u,uN,vN,EX1] = GEN(XS,stdevu,stdevex,n,p,col)
b = ones(col,1);
uN = normrnd(0,1,[1 n]);
vN = normrnd(0,1,[p n]);
u = stdevu*uN;
%EX1 = normrnd(0,stdevex,n,p);
EX1 = stdevex*vN;
EX = [zeros(n,1) EX1'];
covEX = cov(EX1);
X = XS + EX;
y = (XS*b)+ u';
X_seq=ones(n,1);
for j = 1: col
C=X(:,j) - mean(X(:,j),1);
X_seq=[X_seq C];
end
X_seq(:,1:2)=[];
W = X_seq;
z = y -mean(y);
lam = (stdevex^2)/(stdevu^2);
end .........................................................................................................
clc;
rng(156);
global stdevu stdevex Xdiag V vN s rcovar XS1 y n u uN W z lam ;
stdevu = 0.1;
stdevex = 0.3;
n = 10;
p = 1; %%%%%
col = p+1; %%%%%
it = 1;
tmsedls = 0;
XS1 = normrnd(2,1,n,p);%%%%%
XS = [ones(n,1) XS1];%%%%%
rcovar = diag(stdevex*ones(n*p,1));%%%%%%%%
s = ones(n*p,1);
msedls = 0;
[EX,X,y,W,z,lam,u,uN,vN,EX1] = GEN(XS,stdevu,stdevex,n,p,col);
Xdiag = diag(X(:,2));%%%%%
%V = diag(EX1(:,1));%%%%%%%%%%%%
V = diag(vN);%%%%%%%%%%%%
options = optimoptions(@fminunc,'Algorithm','quasi-newton');
b= [1;1];%beta initial value with b0
f = -(((stdevu*y'*uN'+s'*Xdiag*rcovar*V*s)/n)-(y-b(1)*s-XS1*b(2))'*(y-b(1)*s-XS1*b(2))/(2*n)- (s'*V*rcovar^2*V*s)/(2*n))
[b fval] = fminunc(@(b)-(((stdevu*y'*uN'+s'*Xdiag*rcovar*V*s)/n)-(y-b(1)*s-XS1*b(2))'*(y-b(1)*s- XS1*b(2))/(2*n)-(s'*V*rcovar^2*V*s)/(2*n)),b,options);

Antworten (0)

Kategorien

Mehr zu Get Started with Curve Fitting Toolbox finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by