title: How to calculate the standard error of the differences of means in R?
body: <p>The data below shows how vowel durations (vdur) vary by final consonant (finalC). 
  A total of 15 participants produced one item two times. Using <code>ddply</code> function, I got the by-subject means of vowel duration for each item.  </p>
  
  <pre><code>vdur = read.delim("Ldurations.txt")
  Lvdur = vdur[vdur$environment == "A",]
  bySubjmean &lt;- ddply(Lvdur, c("task", "environment", "subject", "item", "finalC", "neutralizedC"), summarize, "bySubj" = mean(vdur), "SD"=sd(vdur))
  
  
  &gt; bySubjmean
  task environment subject   item finalC neutralizedC bySubj          SD
  1      A           A     T01     ip      p            p  90.00  5.37401154
  2      A           A     T01    iph     ph            p  66.05  9.26309883
  3      A           A     T01    kap      p            p  73.80  0.98994949
  4      A           A     T01   kaph     ph            p  48.10 15.98061325
  5      A           A     T01   keot      t            t  79.85  1.20208153
  6      A           A     T01  keoth     th            t  50.65  6.43467171
  7      A           A     T01    mit      t            t  87.45  5.44472222
  8      A           A     T01   mith     th            t  69.55  0.91923882
  9      A           A     T01    muk      k            k  60.25  0.21213203
  10     A           A     T01   mukk     kk            k  42.65  6.01040764 
  11     A           A     T01    nak      k            k 122.10  4.80832611
  12     A           A     T01   nakk     kk            k  66.00  0.98994949
  </code></pre>
  
  <p>I'd like to calculate the standard error of mean vdur differences between finalC "p" and "ph" (and between "t" and "th" and "k" and "kk", respectively). </p>
  
  <p>Although I could easily calculate the mean difference between p and ph (68.42250 - 50.04083 =18.42...) and its SE using ddply(), I was not able to figure out how to calcualte the SE of this mean difference using R codes. </p>
  
  <pre><code>&gt; grandmean &lt;- ddply(bySubjmean, c("task", "environment", "finalC", "neutralizedC"),
     summarize, "grandMean" = mean(bySubj),"SE"=sd(bySubj)/sqrt(length(bySubj)))
  &gt; grandmean
     task environment finalC neutralizedC grandMean       SE
  1     A           A      k            k  84.34133 2.844009
  2     A           A     kk            k  56.51800 1.979809
  3     A           A      p            p  68.42250 2.751806
  4     A           A     ph            p  50.04083 2.208335
  5     A           A      t            t  72.06778 1.975539
  6     A           A     th            t  53.19889 1.819296
  </code></pre>