作者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)