作者buski (......)
看板Statistics
標題[問題] 請問這組R code如果要搬到SAS下該怎麼寫？
時間Fri Nov 17 02:17:56 2006
我是要寫一個跟cross correlation map 有關的語法，我同學用R寫了一組給我
如果我想在SAS下運用，該如何寫呢？ 因為我接下來還要跑Poisson Regression
習慣在SAS下操作
這是某種昆蟲數量跟雨量的關係
VEXF是昆蟲數量
PREC是雨量
precdata=read.csv("D:/cross/prec.csv")
mosdata=read.csv("D:/cross/mos.csv")
vexf=mosdata[,"VEXF"]
prec=precdata[,"PREC"]
mdate=as.Date(mosdata[,"REALDATE"],"%Y/%m/%d")
pdate=as.Date(paste(precdata[,"REALDATE"]%/%10000,
"-",(precdata[,"REALDATE"]%%10000)%/%100,
"-",(precdata[,"REALDATE"]%%100),sep=""))
plot(y=sqrt(vexf),x=mdate)
plot(y=prec,x=pdate)
cor.map=function(x,x.date,z,z.date,Lag1max,Lag2max){
Lag1=0:Lag1max; Nlag1=length(Lag1)
Lag2=0:Lag2max; Nlag2=length(Lag2)
Nd=length(x)
map=matrix(NA,Nlag1,Nlag2)
for(nlag1 in 1:Nlag1){
lag1=Lag1[nlag1]
for(nlag2 in nlag1:Nlag2){
lag2=Lag2[nlag2]
z.avg=rep(NA,Nd)
for(nd in 1:Nd){
z.avg[nd]=mean(z[(z.date>=x.date[nd]-lag2)&(z.date<=x.date[nd]-
lag1)])
}
map[nlag1,nlag2]=cor(x,z.avg)
}
}
list(map=map,lag1=Lag1,lag2=Lag2)
}
cmap=cor.map(x=sqrt(vexf),x.d=mdate,z=prec,z.d=pdate,Lag1=10,Lag2=20)
require(fields)
image.plot(x=cmap$lag1,y=cmap$lag2,xlab="Lag 1",ylab="Lag 2",
z=cmap$map,col=heat.colors(50))
謝謝囉！
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 141.211.10.132
※ 編輯: buski 來自: 141.211.10.132 (11/17 02:19)