用sas的nlp过程求两阶段的负二项分布的最大释然估计

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

用sas的nlp过程求两阶段的负二项分布的最大释然估计。数据是模拟出来的。
程序如下：
data part1(type=est);                                                                                                                                                                                                                                          
      keep _type_  expb0 expb1 theta;                                                                                                                                                                                                                                 
      _type_='parms'; expb0 = 0 ; expb1=1 ; theta=0.5; output;                                                                                                                                                                                                        
      _type_='lb'; expb0 = 1.0e-6 ; expb1 =1.0e-6 ;theta = 1.0e-6;  output;                                                                                                                                                                                          
run;                                                                                                                                                                                                                                                           
proc nlp data=sasuser.com tech=tr inest=part1 outest=opart1    /*b0=0.05,b1=2,gamma=0.2,lamda=5,n=100*/                                                                                                                                                      
outmodel=model1 cov=2 vardef=n pcov phes;                                                                                                                                                                                                                       
max logf;                                                                                                                                                                                                                                                      
parms expb0 expb1 theta;                                                                                                                                                                                                                                             
profile expb0 expb1 theta / alpha = 0.05;                                                                                                                                                                                                                             
y_th = y; x_th = x;                                                                                                                                                                                                                                             
lamda = 5.3600 ;   /* lamda is the mle of the second stage model of possion distribution */                                                                                                                                                                     
b=log(gamma(y_th+theta))-log(gamma(theta))+theta*log(theta) ;                                                                                                
if x_th=1 then c=lamda*expb0*expb1;
if x_th=0 then     c=lamda*expb0;                                                                                                                                                                                                                                          
a=y*log(c)-(y+theta)*log(theta+c);                                                                                                                                                                                                            
s=a+b;                                                                                                                                                                                                                                                         
logf = s;                                                                                                                                                                                                                                                      
run;   
最后算出来的Theta不对呢？请各位高手指点迷津！！
数据产生程序如下（R软件）：
simudata<-function(b0,b1,gamma,lamda,n){
t<-rgamma(n,1/gamma,1,gamma)
m<-rpois(n,lamda)
x<-rbinom(n,1,0.5)
u<-c(1:2)
y<-c(1:2)
p<-c(1:2)
for(i in 1:n)
{ b<-x
if(b<=0)
u<-b0
else
u<-b0*b1
}
for(i in 1:n){
p<-t*u
o<-m
b<-rbinom(1,o,p)
y<-b
}
data<-data.frame(Y=y,M=m,X=x)
print(data)
}
simudata(0.05,2,0.2,5,100)

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：负二项分布 二项分布 sas的 NLP 两阶段 SAS 二项分布 NLP

怎么没人回复呢？

呵呵，我知道哪儿错了~~~~~
参数的参数次幂用npl过程不能求出来，可以用泰勒展开近似求解

cheny600 发表于 2010-5-14 12:02
呵呵，我知道哪儿错了~~~~~
参数的参数次幂用npl过程不能求出来，可以用泰勒展开近似求解

这个也不对。

最后的解决方案是：
用sas的genmod过程提供了负二项分布
最后计算出的结果需要在计算一下。哎，自问自答了，都没人关注。

当然对一些极端的数据，sas的genmod过程提供的算法是不收敛的，如果遇到那种数据就只有自己编程了，建议用s-plus或者R

您好，我想咨询一个问题，我使用9.1版本sas运行nlp时出错，提示“ERROR: 过程 NLP 没有找到”，这是什么原因呢？求帮助！