-
Notifications
You must be signed in to change notification settings - Fork 0
/
cdcMLEtest.m
47 lines (42 loc) · 1.3 KB
/
cdcMLEtest.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
function f=cdcMLEtest(ydataNX,thetaIn)
xdata=1:36;
lx=length(xdata);
mu=[393.0256 999.3054 702.6867 406.0680 177.1958];
sig=[129.6984 329.7708 231.8866 134.0024 58.4746];
nbar=length(mu);
muvec=reshape(repmat(mu,lx,1),nbar*lx,1);
sigvec=reshape(repmat(sig,lx,1),nbar*lx,1);
nbar=size(ydataNX,2);
ydataNX=ydataNX(xdata,:);
hosp=reshape(ydataNX,nbar*lx,1);
thresh=100;
x=zeros(lx*nbar,1);
x(hosp<thresh)=1;
options=optimset('MaxFunEvals',100000);%,'TolX',10^-30,'TolFun',10^-30);
%[theta,~]=mle(hosp,'logpdf',@(hosp,params)lhoodFct1(hosp,params,xdata,lx,nbar,x,muvec,sigvec),'start',thetaIn);
theta=fminsearch(@(params)lhoodFct1(hosp,params,xdata,lx,nbar,x,muvec,sigvec),thetaIn,options);
f=theta;
end
function f=lhoodFct1(hosp,params,xdata,lx,nbar,x,muvec,sigvec)
%{
xdata=1:36;
lx=length(xdata);
mu=[393.0256 999.3054 702.6867 406.0680 177.1958];
sig=[129.6984 329.7708 231.8866 134.0024 58.4746];
nbar=length(mu);
muvec=reshape(repmat(mu,lx,1),nbar*lx,1);
sigvec=reshape(repmat(sig,lx,1),nbar*lx,1);
nbar=size(ydataNX,2);
thresh=100;
x=zeros(lx*nbar,1);
x(hosp<thresh)=1;
%}
ysim=cdcPandemicSimulationW5(params,xdata,0,0,0,243);
y=reshape(ysim,nbar*lx,1);
y(hosp==0)=0;
y=y./hosp;
y(x==1)=NaN;
%muvec(x==1)=NaN;
%sigvec(x==1)=NaN;
f=-nansum(log(normpdf(y,muvec,sigvec)));
end