作者supa666w (灰色哲学)
看板Statistics
标题[程式] Hansen2000门槛模型与R
时间Wed May 11 01:18:03 2011
------------------------------------------------------------------------
[软体程式类别]:R
[程式问题]:程式码
[软体熟悉度]:
新手(不到1个月)
[问题叙述]:
小弟利用R来跑Hansen2000年发表的论文里面的资料
data与code都是作者公布的 他在code里面有设定
可在同质变数或是异质变数的假设下跑回归
code分别是h <-1(假设异质) 或h <-0(假设同质)
在h<-1的情况下我可以跑出完整的报表
可是当我输入h <-0再输入回归式时就出现
错误在b1 - ser1 * z : 非调和阵列
可是我只是在h那里改0或1而已 其他都设定好了也没变动
不知道哪位前辈可以为小弟解惑?
第一次在板上po文 无理之处烦请见谅
[程式范例]:
# Load procedures and data #
source("C:/r/thr_het.R")
source("c:/r/thr_est.R")
data <- read.table("c:/r/dur_john.dat")
# Exclude missing variables and oil states #
k <- ncol(data)
indx <- as.matrix(1-(data[,5]== -999))%*%matrix(c(1),1,k)
data <- as.matrix(data[indx>0])
data <- matrix(data,nrow=nrow(data)/k,ncol=k)
indx <- as.matrix(1-(data[,6]== -999))%*%matrix(c(1),1,k)
data <- as.matrix(data[indx>0])
data <- matrix(data,nrow=nrow(data)/k,ncol=k)
indx <- as.matrix(1-(data[,10]== -999))%*%matrix(c(1),1,k)
data <- as.matrix(data[indx>0])
data <- matrix(data,nrow=nrow(data)/k,ncol=k)
indx <- as.matrix(1-(data[,11]== -999))%*%matrix(c(1),1,k)
data <- as.matrix(data[indx>0])
data <- matrix(data,nrow=nrow(data)/k,ncol=k)
indx <- as.matrix(data[,2]== 1)%*%matrix(c(1),1,k)
data <- as.matrix(data[indx>0])
data <- matrix(data,nrow=nrow(data)/k,ncol=k)
# Make data transformations #
diff <- log(data[,6])-log(data[,5])
q <- data[,5]
gdp60 <- log(q)
iony <- log(data[,9]/100)
pgro <- log(data[,8]/100+.05)
sch <- log((data[,10])/100)
lit <- data[,11]
dat <- cbind(diff,gdp60,iony,pgro,sch,q,lit)
# Program switches #
rep <- 1000
h <- 1 ------------->问题所在
na <- rbind("GNP_Gwth","GDP_1960","Inv/GDP","Pop_Gwth","School","GDP_1960",
"Literacy")
dum <- rbind(2,3,4,5)
# Estimate First Sample Split, Using Output as Threshold #
qhat1 <- thr_est(dat,na,1,dum,6,h)
---------------------------------------------------------------------------
--
※ 发信站: 批踢踢实业坊(ptt.cc)
◆ From: 123.194.244.102
1F：→ gsuper:as.matrix(1-(data[,5]== -999))%*%matrix(c(1),1,k) 05/11 01:30
2F：→ gsuper:这行非常诡异 05/11 01:30
3F：→ gsuper:首先 data[,5] 是向量 , data[,5]== -999 是0,1的向量 05/11 01:31
4F：→ gsuper:1-(data[,5]== -999) 变成 0,-1 的向量 05/11 01:31
5F：→ gsuper:再来是 matrix(c(1),1,k) 是 k个1 05/11 01:32
6F：→ gsuper:所以在 0,-1的向量中 , 找 1 这个数值,sum(indx)恒等於0 05/11 01:33
7F：→ gsuper:我想你会不会是程式码哪里贴错了? 05/11 01:33
8F：→ gsuper:或是非调和阵列的原因 , 出自於 data 中有missing variable 05/11 01:35
9F：→ gsuper:因为第二段程式真的很怪 05/11 01:35
10F：→ supa666w:第二段程式码的标题是去除Missing variables 05/11 01:40
11F：→ supa666w:所以是没问题的 我贴的code是作者公布的 网址如下 05/11 01:41
12F：推 gsuper:我发现我把 %*% 看成 %in% , 我再想想看好了 05/11 01:44
14F：→ supa666w:谢谢你热心的回覆:) 05/11 01:44
15F：推 Wush978:看起来问题在 thr_est.r 里面有一段 05/11 10:49
16F：→ Wush978: if (h==0){ 05/11 10:49
17F：→ Wush978:ser1 <- as.matrix(sqrt(diag(mi1)*(t(e1)%*%e1)/(nrow(y1 05/11 10:49
18F：→ Wush978:厄... 第239行 05/11 10:49
19F：→ Wush978:(nrow(y1)-k) 05/11 10:49
20F：→ Wush978:这时候y1 的资料型态已经不再是matrix 05/11 10:50
21F：→ Wush978:nrow(y1) return numeric(0) 05/11 10:50
22F：→ Wush978:逻辑的部份我不细看了, 请自行解决吧 :) 程式的问题在这 05/11 10:50
23F：→ supa666w:谢谢你的回覆 我会试试看的! 05/12 10:49