为什么不早出现这个帖子啊，555，不知道parse这个函数，期中考试题都没做好[s:15]
感觉求mle的时候这个函数非常好使，看了这个帖子又做了下考试题，不知道大家有没有做过类似的练习，多多指教[s:13]
<br />
#Newton<br />
mle<-function(x,theta0,epsilon=1.0e-5,iterations=100){<br />
n=length(x)<br />
e=1<br />
d=epsilon<br />
N=0<br />
partial="(1+exp(theta-x[1]))^-1"<br />
for(i in 2:n)<br />
partial=paste(partial,"+(1+exp(theta-x[",i,"]))^-1",sep="")<br />
dltheta=paste("2*","(",partial,")",-n,sep="")<br />
f<-function(theta){eval(parse(text=dltheta))}</p>
<p> while(e>epsilon){<br />
N=N+1<br />
if(N>iterations) cat("not converge after",N,"iterations")<br />
theta1=theta0-f(theta0)*d/(f(theta0+d)-f(theta0))<br />
e=abs(theta1-theta0)<br />
theta0=theta1<br />
}<br />
return(list(theta=theta0,n=N))<br />
}<br />
#================= test =======================<br />
set.seed(555)<br />
x=rlogis(100,3,1)<br />
mle(x,0)<br />
##result:<br />
#$theta<br />
#[1] 2.878581<br />
#$n<br />
#[1] 5<br />
#================================================================<br />
#Bisection<br />
mle<-function(x,theta0,theta1,epsilon=1.0e-5){<br />
n=length(x)<br />
e=1<br />
partial="(1+exp(theta-x[1]))^-1"<br />
for(i in 2:n)<br />
partial=paste(partial,"+(1+exp(theta-x[",i,"]))^-1",sep="")<br />
dltheta=paste("2*","(",partial,")",-n,sep="")<br />
f<-function(theta){eval(parse(text=dltheta))}</p>
<p> y0=f(theta0)<br />
y1=f(theta1)<br />
if(y0*y1>0) return ("wrong initial values")<br />
error=abs(theta1-theta0)<br />
N=0<br />
while(error>epsilon){<br />
N=N+1<br />
error=error/2<br />
theta2=(theta0+theta1)/2<br />
y2=f(theta2)<br />
if(y0*y2<0){<br />
theta1=theta2; y1=y2}<br />
else {<br />
theta0=theta2; y0=y2}<br />
}<br />
return(list(theta=theta2, n=N))<br />
}</p>
<p>#================ test ========================<br />
set.seed(555)<br />
x=rlogis(100,3,1)<br />
mle(x,1,5)<br />
##result<br />
#$theta<br />
#[1] 2.878578<br />
#$n<br />
#[1] 19</p>
<p>
</p>
[attachment=210342,629]
