Rコーディング
ひたすら、Rのコードを打ち込む一日。
とりあえず、今日書いたコードをさらしあげ。
forの多重ネストばっかで、すごくいけてない（つまり、すんごく遅い）けど、いろいろ考える気力がおきましぇん。
明日くらいからbugsコード書き始めたい。
[sourcecode language=’cpp’]library(R2WinBUGS)
library(coda)
library(VGAM)
setwd(“d:/【研究】強度間伐関連/強度間伐データ”)
first.dat <- read.csv("１期目実生調査.csv",header=T)
first.dat$No <- as.character(first.dat$No)ed
first.dat$Nokabu <- as.character(first.dat$Nokabu)
first.dat$sp <- as.character(first.dat$sp)
first.dat$nos <- as.numeric(first.dat$nos)
first.sp <- unique(first.dat$sp)
first.sp <- first.sp[first.sp !="nothing"]
first.name.dat <- sort(as.character(unique(first.dat$name)))
second.dat<-read.csv("２期目実生調査.csv",header=T)
second.dat$No <- as.character(second.dat$No)
second.dat$Nokabu <- as.character(second.dat$Nokabu)
second.dat$sp <- as.character(second.dat$sp)
second.dat$nos <- as.numeric(second.dat$nos)
second.sp <- unique(second.dat$sp)
second.sp <- second.sp[second.sp !="nothing"]
second.name.dat <- sort(as.character(unique(second.dat$name)))
sp <- append(first.sp,second.sp)
sp <- (sort(unique(sp)))
##毎木データをプロットごとのデータフレームに分ける
##１期目のデータのループ
#############################################################
##各プロットの毎木データをそれぞれのデータフレームに分ける。
x <- 0
subplot.name.dat <- NULL
for( i in 1:length(first.name.dat)){
x <- x+1
for( subID in 1:4){
each.dataframe <- subset(first.dat,first.dat$name == first.name.dat[x] & first.dat$subplot == subID )
##assign()でプロットＩＤをオブジェクト名に当てる。
each.name <- paste(first.name.dat[x],subID,sep="")
subplot.name.dat <- append(subplot.name.dat,each.name)
assign(each.name,each.dataframe)
}
}
##種・プロットのクロス表を作成
##１期目のデータのループ
############################
##空のデータフレーム作る
first.sp.dat <- NULL
for(mat in 1:length(subplot.name.dat)-1){
first.sp.dat <- rbind(first.sp.dat,c(rep(0,length=length(sp))))
}
rownames(first.sp.dat) <- subplot.name.dat
colnames(first.sp.dat) <- sp
##各種をカウントしていくループ
for( k in 1:length(subplot.name.dat)){
mod.dat <- NULL
for( m in 1:length(get(subplot.name.dat[k])$name)){
mod.dat <- rbind(mod.dat,get(subplot.name.dat[k])[m,])
}
for(n in 1:length(sp)){
for(p in 1:length(mod.dat$name)){
if(sp[n] == mod.dat$sp[p]){
first.sp.dat[k,n] <- first.sp.dat[k,n] + mod.dat$nos[p]
}else{
}
}
}
}
##毎木データをプロットごとのデータフレームに分ける
##２期目のデータのループ
#############################################################
##各プロットの毎木データをそれぞれのデータフレームに分ける。
x <- 0
subplot.name.dat <- NULL
for( i in 1:length(name.dat)){
x <- x+1
for( subID in 1:4){
each.dataframe <- subset(second.dat,second.dat$name == name.dat[x] & second.dat$subplot == subID )
##assign()でプロットＩＤをオブジェクト名に当てる。
each.name <- paste(second.name.dat[x],subID,sep="")
subplot.name.dat <- append(subplot.name.dat,each.name)
assign(each.name,each.dataframe)
}
}
##種・プロットのクロス表を作成
##２期目のデータのループ
############################
##空のデータフレーム作る
second.sp.dat <- NULL
for(mat in 1:length(subplot.name.dat)-1){
second.sp.dat <- rbind(second.sp.dat,c(rep(0,length=length(sp))))
}
rownames(second.sp.dat) <- subplot.name.dat
colnames(second.sp.dat) <- sp
##各種をカウントしていくループ
for( k in 1:length(subplot.name.dat)){
mod.dat <- NULL
for( m in 1:length(get(subplot.name.dat[k])$name)){
if(get(subplot.name.dat[k])$Nokabu[m]==get(subplot.name.dat[k])$No[m] | is.na(get(subplot.name.dat[k])$No[m])){
mod.dat <- rbind(mod.dat,get(subplot.name.dat[k])[m,])
}else{
}
}
for(n in 1:length(sp)){
for(p in 1:length(mod.dat$name)){
if(sp[n] == mod.dat$sp[p]){
second.sp.dat[k,n] <- second.sp.dat[k,n] + mod.dat$nos[p]
}else{
}
}
}
}
diff.sp.dat <- second.sp.dat - first.sp.dat
##種の出現頻度を算出
col.count <- NULL
threshold <- 3
for(freq.sp in 1:length(diff.sp.dat[1,])){
if(length(diff.sp.dat[,freq.sp][diff.sp.dat[,freq.sp] != 0]) >threshold ){
col.count <- append(col.count,freq.sp)
}else{
}
}
mod.sp.dat <- diff.sp.dat[,col.count][/sourcecode]
-
前の記事
最近の独り言 on 2009-03-15 2009.03.15
-
次の記事
MCMCの神様が降りてきた 2009.03.17