2

複数の応答変数と 3 つの処理を含むデータセットがあります。治療 2 は治療 1 内にネストされ、治療 3 は治療 2 内にネストされます。簡単にするために、3 つの応答変数のみを示しています。これを 22 以上の応答変数で実行したいと思います。そのうちの 3 つがデモ テーブルに表示されています。

私の目的:

  1. 治療の組み合わせに基づいて応答変数がどのように変化するかを視覚化します。1 つの応答変数でこれを実行するスクリプトを作成しました。このコードをコピーして他の列を実行していますが、これは非常に粗雑な方法です。これが私の 2 番目の目的につながります。
  2. 次のスクリプトを自動化または変更して、列を自動的にループし、目的の表とグラフを生成できるようにします。

デモデータ: demo.table

これが私のスクリプトです：

library(doBy)
length2 <- function (x, na.rm=FALSE) {
    if (na.rm) sum(!is.na(x))
    else       length(x)
}
 attach (demo)
cdataNA <- summaryBy(tyr ~ spp + wat + ins, data=demo, FUN=c(length2,mean,sd), na.rm=TRUE)
# Rename column change.length to just N
names(cdataNA)[names(cdataNA)=="tyr.length2"] <- "N"
# Calculate standard error of the mean
cdataNA$tyr.SE <- cdataNA$tyr.sd / sqrt(cdataNA$N)
cdataNA
# Now create a barplot using ggplot2
library(ggplot2)
a <- ggplot(cdataNA, aes(x = wat, y = tyr.mean, fill = ins))
b <- a + geom_bar(stat = "identity", position = "dodge") + facet_grid (~ spp)
# Now put errorbars.
c <- b + geom_errorbar(aes(ymin=tyr.mean-tyr.SE, ymax=tyr.mean+tyr.SE), 
                  width=.2,                    # Width of the error bars
                  position=position_dodge(.9)) + 
xlab ("wat") + 
ylab ("tyr (PA/PA std)")
c

## esc
library(doBy)
length2 <- function (x, na.rm=FALSE) {
    if (na.rm) sum(!is.na(x))
    else       length(x)
}
cdataNA1 <- summaryBy(esc ~ spp + wat + ins, data=demo, FUN=c(length2,mean,sd), na.rm=TRUE)
# Rename column change.length to just N
names(cdataNA1)[names(cdataNA1)=="esc.length2"] <- "N"
# Calculate standard error of the mean
cdataNA1$esc.SE <- cdataNA1$esc.sd / sqrt(cdataNA1$N)
cdataNA1
# Now create a barplot using ggplot2
library(ggplot2)
a1 <- ggplot(cdataNA1, aes(x = wat, y = esc.mean, fill = ins))
b1 <- a1 + geom_bar(stat = "identity", position = "dodge") + facet_grid (~ spp)
# Now put errorbars.
c1 <- b1 + geom_errorbar(aes(ymin=esc.mean-esc.SE, ymax=esc.mean+esc.SE), 
                  width=.2,                    # Width of the error bars
                  position=position_dodge(.9)) + 
xlab ("wat") + 
ylab ("esc (PA/PA std)")
c1

tyr の結果テーブル:

  spp  wat ins N tyr.mean      tyr.sd      tyr.SE
1  Bl High  No 4 0.305325 0.034102041 0.017051020
2  Bl High Yes 5 0.186140 0.045165894 0.020198802
3  Bl  Low  No 5 0.310540 0.061810096 0.027642315
4  Bl  Low Yes 5 0.202840 0.029034944 0.012984822
5 Man High  No 4 0.122725 0.075867005 0.037933503
6 Man High Yes 5 0.081800 0.013463469 0.006021046
7 Man  Low  No 5 0.079880 0.009569587 0.004279650
8 Man  Low Yes 4 0.083550 0.018431947 0.009215973

esc の結果グラフ: esc のデモ図

したがって、すべてが機能しますが、ワークフローを妨げるかなりの手作業が必要です。自動化を達成することは素晴らしいことです。

前もって感謝します。

4

1 に答える 1

0

次の 2 行だけでデータを整理できます。

melt.dta <- melt(dta, id.vars = c("spp", "wat", "ins"), measure.vars = "tyr")
cast(melt.dta, spp + wat + ins ~ ., 
     function (x) c("N" = sum(!is.na(x)), 
                    "mean" = mean(x, na.rm = TRUE),                   
                    "sd" = sd(x, na.rm = TRUE),
                    "se" = sd(x, na.rm = TRUE)/sqrt(sum(!is.na(x)))))

戻り値:

  spp  wat ins N   mean      sd      se
1  Bl High  No 4 0.3053 0.03410 0.01705
2  Bl High Yes 5 0.1861 0.04517 0.02020
3  Bl  Low  No 5 0.3105 0.06181 0.02764
4  Bl  Low Yes 5 0.2028 0.02903 0.01298
5 Man High  No 4 0.1227 0.07587 0.03793
6 Man High Yes 5 0.0818 0.01346 0.00602
7 Man  Low  No 5 0.0799 0.00957 0.00428
8 Man  Low Yes 4 0.0835 0.01843 0.00922
于 2012-10-20T02:20:03.233 に答える