続いて、本来の目的である個人ベースの言及数の解析に入ります。こんな感じのスクリプトで解析してみることにします。
data.stat <- function(d, pickup=seq_len(ncol(d))){ d <- d[,pickup] lm.ba <- apply(d, 1, function(v) lm(as.numeric(v) ~ pickup)[["coefficients"]]) data.frame(name=rownames(d), total=apply(d, 1, sum), a=lm.ba[2,], b=lm.ba[1,], d, row.names=seq_len(nrow(d))) } length.header <- 4 plot.countsy <- function(x, pickup=seq_len(length(x)-length.header), vvalue){ intercept <- x[["b"]] slope <- x[["a"]] count.data <- as.numeric(x[-(1:length.header)]) label <- substr(names(x)[-(1:length.header)], 1, 3) range.y <- c(pickup*slope + intercept, count.data) # plot data plot(pickup, count.data, xaxt="n", ylim=c(min(range.y), max(range.y)), type="o", pch=16, col=palette()[1], lwd=2, ann=F) # plot linear model lines(pickup, pickup*slope+intercept, lwd=2, col=palette()[4]) # write labels axis(1, pickup, label) if(!missing(vvalue)) abline(v=vvalue, col="#CCCCCC", lty=2) }
前回の関数と合わせて、こんな感じ。
d <- read.sy.data("2008_sy.txt") d1 <- d[,rev(seq(ncol(d), by=-1, length=20))] # exclude 6,7,8,12(2007) and 1,6(2008) pickup <- seq_len(ncol(d1))[-c(2:4,8,9,14)] x <- d2[(d2$name == "%B9%E2%B3%C0%BA%CC%CD%DB"),] plot.countsy(x, pickup, 9)
・・・とすると、こんなグラフが出てきます。
data.stat()は、言及数の合計値と線形回帰の傾きと切片を出している関数。この場合は2007年も含めて計算しているので、2008年だけ取りたければ適当にフィルタリングする必要があります。
で、例えば合計値でソートしたいときは
> head(d2[sort.list(d2$total, dec=T), 1:2]) name total 3 %A4%A2%A4%EA%A4%B9 19596.804 775 %CA%BF%CC%EE%B0%BD 7830.758 476 %BF%E5%BC%F9%C6%E0%A1%B9 7445.707 675 %C5%C4%C2%BC%A4%E6%A4%AB%A4%EA 6814.367 637 %C3%E6%C0%EE%E6%C6%BB%D2 5820.329 797 %CB%D9%B9%BE%CD%B3%B0%E1 5171.893
こんな感じ。
やっぱり
Rはお手軽すぎる。同じことをPerlでやろうと思ったらどれだけめんどくさいか・・・。