# (Intercept)      0.0929170  0.0050641  18.348   <2e-16 ***
# logAllBoundFreq -0.0015581  0.0006542  -2.382   0.0175 *
# ---
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1   1
#
# Residual standard error: 0.03927 on 645 degrees of freedom
# Multiple R-squared:  0.008718,  Adjusted R-squared:  0.007181
# F-statistic: 5.672 on 1 and 645 DF,  p-value: 0.01752
c4.lm <- lm (lenMorph ~ logBaseFreq + morph, data = d)
summary(c4.lm)
# Call:
# lm(formula = lenMorph ~ logBaseFreq + morph, data = d)
#
# Residuals:
#       Min        1Q    Median        3Q       Max
# -0.073569 -0.022174 -0.006519  0.013105  0.153842
#
# Coefficients:
#               Estimate Std. Error t value Pr(>|t|)
# (Intercept)  0.1182427  0.0052940  22.335  < 2e-16 ***
# logBaseFreq -0.0018144  0.0006623  -2.740  0.00632 **
# morphplural -0.0181679  0.0045790  -3.968 8.08e-05 ***
# morph3rdsg  -0.0313100  0.0047518  -6.589 9.26e-11 ***
# morphGEN    -0.0302321  0.0047188  -6.407 2.89e-10 ***
# morphhas    -0.0427487  0.0060361  -7.082 3.76e-12 ***
# morphis     -0.0345364  0.0047337  -7.296 8.83e-13 ***
# morphPL-GEN -0.0366526  0.0078438  -4.673 3.63e-06 ***
# ---
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1   1
#
# Residual standard error: 0.0354 on 639 degrees of freedom
# Multiple R-squared:  0.2023,    Adjusted R-squared:  0.1935
# F-statistic: 23.15 on 7 and 639 DF,  p-value: < 2.2e-16
anova(c1.lm, c4.lm)
# Analysis of Variance Table
#
# Model 1: lenMorph ~ morph
# Model 2: lenMorph ~ logBaseFreq + morph
#   Res.Df     RSS Df Sum of Sq      F   Pr(>F)
# 1    640 0.81007
# 2    639 0.80067  1 0.0094043 7.5054 0.006323 **
# ---
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1   1
c5.lm <- lm (lenMorph ~ morph + logAllBoundFreq , data = d)
summary(c5.lm)
#  Call:
# lm(formula = lenMorph ~ morph + logAllBoundFreq, data = d)
#
# Residuals:
#      Min       1Q   Median       3Q      Max
# -0.07380 -0.02217 -0.00679  0.01355  0.15442
#
# Coefficients:
#                   Estimate Std. Error t value Pr(>|t|)
# (Intercept)      0.1166688  0.0050307  23.191  < 2e-16 ***
# morphplural     -0.0191062  0.0045148  -4.232 2.66e-05 ***
# morph3rdsg      -0.0352665  0.0043737  -8.063 3.67e-15 ***
# morphGEN        -0.0332962  0.0045397  -7.335 6.77e-13 ***
# morphhas        -0.0476928  0.0057451  -8.301 6.13e-16 ***
# morphis         -0.0391756  0.0044174  -8.868  < 2e-16 ***
# morphPL-GEN     -0.0424526  0.0079145  -5.364 1.14e-07 ***
# logAllBoundFreq -0.0015908  0.0006193  -2.569   0.0104 *
# ---
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1   1
#
# Residual standard error: 0.03542 on 639 degrees of freedom
# Multiple R-squared:  0.2011,    Adjusted R-squared:  0.1924
# F-statistic: 22.98 on 7 and 639 DF,  p-value: < 2.2e-16
anova(c1.lm, c5.lm)
# Analysis of Variance Table
#
# Model 1: lenMorph ~ morph
# Model 2: lenMorph ~ morph + logAllBoundFreq
#   Res.Df     RSS Df Sum of Sq      F  Pr(>F)
# 1    640 0.81007
# 2    639 0.80180  1 0.0082782 6.5974 0.01044 *
# ---
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1   1
# we improve the morph model if we add baseFreq, but not if we add AllBoundFreq
c4.lmInt <- lm (lenMorph ~ morph * logAllBoundFreq , data = d)
summary(c4.lmInt)
# ns
visreg(c4.lmInt, "logAllBoundFreq", by="morph")
c5.lmInt <- lm (lenMorph ~ morph * logBaseFreq , data = d)
summary(c5.lmInt)
# ns
visreg(c5.lmInt, "logBaseFreq", by="morph")
summary(d$logTrigramFreq)
#    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's
#   0.000   0.000   0.000   1.127   1.609   8.379     156
summary(d$logLBigramFreq)
#    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's
#   0.000   0.000   2.773   3.111   5.204  10.110      64
summary(d$logRBigramFreq)
#    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's
#   0.000   0.000   1.609   2.540   4.331   9.884      98
summary(d$logTrigramFreq)
#    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's
#   0.000   0.000   0.000   1.127   1.609   8.379     156
par(mfrow=c(1,3))
densityplot(d$logTrigramFreq)
densityplot(d$logLBigramFreq)
densityplot(d$logRBigramFreq)
par(mfrow=c(1,1))
plot(d$logLBigramFreq, d$logRBigramFreq)
plot1 <- ggplot(d, aes(logLBigramFreq,logRBigramFreq))
plot1 + geom_point() + labs(x="log left bigram frequency", y="log right bigram frequency") + geom_smooth()
plot1 + geom_point() + labs(x="log left bigram frequency", y="log right bigram frequency") + geom_smooth(method="lm")
######### some descriptive stats
tapply(d$lenMorph, d$morph, mean)
# s          plural     3rdsg      GEN        has        is         PL-GEN
# 0.10547199 0.08396957 0.06885465 0.07163406 0.05761365 0.06625264 0.06644534
range(tapply(d$lenMorph, d$morph, mean))
# [1] 0.05761365 0.10547199
0.10547199-0.05761365
#[1] 0.04785834
tapply(d$lenMorph, d$isVoiced, mean)
# voiced     unvoiced
# 0.05764802 0.08483373
d$morph <- relevel(d$morph, "plural")
d$morph <- relevel(d$morph, "s")
### now checking collinearity of duration and syllable number
summary(d$sylNum)
#    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
#    1.00    1.00    2.00    1.72    2.00    6.00
summary(d$lenMorph)
#    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
# 0.01936 0.05376 0.07159 0.08143 0.09802 0.23700
cor.test(d$sylNum, d$lenBase)
#         Pearson's product-moment correlation
#
# data:  d$sylNum and d$lenBase
# t = 16.95, df = 645, p-value < 2.2e-16
# alternative hypothesis: true correlation is not equal to 0
# 95 percent confidence interval:
#  0.4994099 0.6062599
# sample estimates:
#       cor
# 0.5551209
syll.lm <- lm (lenMorph ~ morph + sylNum , data = d)
summary(syll.lm)
# Call:
#   lm(formula = lenMorph ~ morph + sylNum, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.076650 -0.022775 -0.007314  0.013946  0.157357
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
# (Intercept)  0.098566   0.003825  25.766  < 2e-16 ***
#   morphplural -0.020808   0.004430  -4.697 3.23e-06 ***
#   morph3rdsg  -0.034602   0.004426  -7.818 2.21e-14 ***
#   morphGEN    -0.034694   0.004552  -7.622 9.06e-14 ***
#   morphhas    -0.048434   0.005754  -8.418 2.52e-16 ***
#   morphis     -0.041206   0.004498  -9.162  < 2e-16 ***
#   morphPL-GEN -0.039326   0.007808  -5.037 6.16e-07 ***
#   sylNum       0.004042   0.001688   2.394   0.0169 *
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03545 on 639 degrees of freedom
# Multiple R-squared:  0.2001,  Adjusted R-squared:  0.1913
# F-statistic: 22.83 on 7 and 639 DF,  p-value: < 2.2e-16
syll.lm1 <- lm (lenMorph ~ morph + actSylNum , data = d)
summary(syll.lm1)
# Call:
#   lm(formula = lenMorph ~ morph + actSylNum, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.077172 -0.022646 -0.007307  0.014104  0.157391
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
#   (Intercept)  0.097071   0.003866  25.106  < 2e-16 ***
#   morphplural -0.020864   0.004418  -4.723 2.86e-06 ***
#   morph3rdsg  -0.034328   0.004410  -7.784 2.85e-14 ***
#   morphGEN    -0.035251   0.004556  -7.737 3.99e-14 ***
#   morphhas    -0.048162   0.005739  -8.392 3.07e-16 ***
#   morphis     -0.041186   0.004466  -9.223  < 2e-16 ***
#   morphPL-GEN -0.039410   0.007793  -5.057 5.57e-07 ***
#   actSylNum    0.005051   0.001769   2.855  0.00444 **
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03538 on 639 degrees of freedom
# Multiple R-squared:  0.2031,  Adjusted R-squared:  0.1943
# F-statistic: 23.26 on 7 and 639 DF,  p-value: < 2.2e-16
pairscor.fnc(d[, c("morph", "lenMorph", "lenBase", "actSylNum", "sylNum")])
# sylNum and lenBase correlate with each other, but only lenMorph correlates with lenBase
syll.lm2 <- lm (lenMorph ~ morph + lenBase, data = d)
summary(syll.lm2)
# Call:
#   lm(formula = lenMorph ~ morph + lenBase, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.098409 -0.021010 -0.006282  0.013546  0.150424
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
#   (Intercept)  0.073641   0.004409  16.704  < 2e-16 ***
#   morphplural -0.019125   0.004213  -4.539 6.74e-06 ***
#   morph3rdsg  -0.026931   0.004283  -6.288 5.97e-10 ***
#   morphGEN    -0.035223   0.004319  -8.156 1.84e-15 ***
#   morphhas    -0.040622   0.005532  -7.343 6.37e-13 ***
#   morphis     -0.039291   0.004204  -9.346  < 2e-16 ***
#   morphPL-GEN -0.037646   0.007426  -5.069 5.23e-07 ***
#   lenBase      0.090045   0.010480   8.592  < 2e-16 ***
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03371 on 639 degrees of freedom
# Multiple R-squared:  0.2765,  Adjusted R-squared:  0.2686
# F-statistic: 34.88 on 7 and 639 DF,  p-value: < 2.2e-16
syll.lm3 <- lm (lenMorph ~ lenBase, data = d)
summary(syll.lm3)
# Call:
#   lm(formula = lenMorph ~ lenBase, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.080639 -0.025374 -0.007832  0.016596  0.161039
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
#   (Intercept) 0.048810   0.003887  12.558   <2e-16 ***
#   lenBase     0.099202   0.010954   9.056   <2e-16 ***
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03716 on 645 degrees of freedom
# Multiple R-squared:  0.1128,  Adjusted R-squared:  0.1114
# F-statistic: 82.02 on 1 and 645 DF,  p-value: < 2.2e-16
syll.lm4 <- lm (lenMorph ~  actSylNum , data = d)
summary(syll.lm4)
#
# Call:
#   lm(formula = lenMorph ~ actSylNum, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.063280 -0.028006 -0.009147  0.016517  0.158087
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
#   (Intercept) 0.075140   0.003492  21.516   <2e-16 ***
#   actSylNum   0.003751   0.001867   2.009    0.045 *
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03932 on 645 degrees of freedom
# Multiple R-squared:  0.006217,  Adjusted R-squared:  0.004677
# F-statistic: 4.035 on 1 and 645 DF,  p-value: 0.04497
syll.lm5 <- lm (lenMorph ~  lenBase + actSylNum , data = d)
summary(syll.lm5)
#
# Call:
#   lm(formula = lenMorph ~ lenBase + actSylNum, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.084078 -0.025293 -0.007189  0.015953  0.159464
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
#   (Intercept)  0.052931   0.003947  13.410  < 2e-16 ***
#   lenBase      0.135784   0.013657   9.942  < 2e-16 ***
#   actSylNum   -0.009631   0.002200  -4.378  1.4e-05 ***
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03664 on 644 degrees of freedom
# Multiple R-squared:  0.1385,  Adjusted R-squared:  0.1358
# F-statistic: 51.75 on 2 and 644 DF,  p-value: < 2.2e-16
# We find a suppression effect!
syll.lm6 <- lm (lenMorph ~  morph + lenBase + actSylNum , data = d)
summary(syll.lm6)
# Call:
#   lm(formula = lenMorph ~ morph + lenBase + actSylNum, data = d)
#
# Residuals:
#   Min        1Q    Median        3Q       Max
# -0.099515 -0.020221 -0.006575  0.013364  0.149114
#
# Coefficients:
#   Estimate Std. Error t value Pr(>|t|)
#   (Intercept)  0.075605   0.004452  16.981  < 2e-16 ***
#   morphplural -0.019284   0.004194  -4.598 5.16e-06 ***
#   morph3rdsg  -0.027241   0.004265  -6.387 3.27e-10 ***
#   morphGEN    -0.034007   0.004324  -7.864 1.59e-14 ***
#   morphhas    -0.038678   0.005556  -6.961 8.41e-12 ***
#   morphis     -0.037184   0.004262  -8.725  < 2e-16 ***
#   morphPL-GEN -0.036924   0.007397  -4.991 7.74e-07 ***
#   lenBase      0.110146   0.012955   8.502  < 2e-16 ***
#   actSylNum   -0.005453   0.002084  -2.617  0.00908 **
#   ---
#   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#
# Residual standard error: 0.03356 on 638 degrees of freedom
# Multiple R-squared:  0.2842,  Adjusted R-squared:  0.2752
# F-statistic: 31.66 on 8 and 638 DF,  p-value: < 2.2e-16
# suppression effect remains
# bottom line:
# we model with lenBase (excluding SylNum) and with baseFreq (excluding allBoundFreq)
####
# Box Cox transformation to correct for skewed distribtion of lenMorph
####
tmp.lm <- lm(lenMorph ~
morph +
position +
consonants +
sylSec +
logBaseFreq +
isVoiced +
baseRep +
lenBase+
logLBigramFreq +
logRBigramFreq
, data = d)
bc.d <- boxcox (tmp.lm)
Exponent <- bc.d$x [which.max (bc.d$y)]
Exponent
# [1] 0.1010101
d$lenMorphBC <- d$lenMorph^Exponent
# we check for strange data points
densityplot(d$sylSec)
# we need to remove two points
d[d$sylSec > 15,]
d <- d[d$sylSec < 15,]
dim(d)
# [1] 645  50
summary(d$morph)
sum(summary(d$morph))
#      s plural  3rdsg    GEN    has     is PL-GEN   NA's
#    196     95    100     88     47     95     23      1
# The NA has to be removed
d[is.na(d$morph==T),]
d1 <- d[is.na(d$morph)==F,]
dim(d1)
# 644 50
d <- d1
t <- table(d$morph)
t
# s plural  3rdsg    GEN    has     is PL-GEN
# 196     95    100     88     47     95     23
sum(t)
# 644
# number of morphemic S
644-196
# 448
table(d$consonants)
# 0   1   2   3
# 325 259  58   2
summary(d$consonants)
#   0   1   2   3
# 325 259  58   2
sum(summary(d$consonants))
# [1] 644
summary(d$isVoiced)
#   voiced unvoiced
#       81      563
summary(d$position)
#  final middle
#     97    547
summary(d$boundary)
#  no yes
# 418 226
summary(d$succSound)
#    $    _    {    1    2    3    5    6   aa    b    d    D    E    f    g    h    H    i    I    j    J    k    l    m    n   ow    p    P    Q    r
#    7    5   11    2   11    7    4   10    1   33   13   20   10   26   24   25   18    2   48   11    5   20   19   16   17    1   11    2    8   14
#    s    S    t    T    u    v    V    w    z NA's
#   38    5   22    4    1    2   49   22    3   97
sum(summary(d$succSound))
# [1] 644
644-97
# [1] 547
unique(d$succSound)
#  [1] p    b    d    f    w    h    g    k    v    s    J    $    l    {    H    r    <NA> 2    3    m    I    j    6    V    n    t    E    5    _
# [30] Q    D    u    z    T    P    i    S    aa   1    ow
# Levels: $ _ { 1 2 3 5 6 aa b d D E f g h H i I j J k l m n ow p P Q r s S t T u v V w z
summary(d$morph)
#      s plural  3rdsg    GEN    has     is PL-GEN
#    196     95    100     88     47     95     23
summary(d$boundary)
#  no yes
# 418 226
names(d)
pairscor.fnc(d[ , c("morph", "position", "sylSec", "logBaseFreq", "isVoiced", "logLBigramFreq", "logRBigramFreq", "lenBase",  "consonants")])
pairscor.fnc(d[ , c("morph", "position", "sylSec", "logBaseFreq", "isVoiced", "lenBase")])
aovmodel <- aov(lenMorph ~ morph,d)
anova(aovmodel)
# Analysis of Variance Table
#
# Response: lenMorph
#            Df  Sum Sq  Mean Sq F value    Pr(>F)
# morph       6 0.18767 0.031279  24.864 < 2.2e-16 ***
# Residuals 637 0.80135 0.001258
# ---
# Signif. codes:  0 *** 0.001 ** 0.01 * 0.05 . 0.1   1
TukeyHSD(aovmodel)
# Tukey multiple comparisons of means
# 95% family-wise confidence level
#
# Fit: aov(formula = lenMorph ~ morph, data = d)
#
# $morph
# diff         lwr           upr     p adj
# plural-s      -0.0210569652 -0.03417186 -7.942069e-03 0.0000517
# 3rdsg-s       -0.0361718792 -0.04906405 -2.327971e-02 0.0000000
# GEN-s         -0.0333924711 -0.04685411 -1.993083e-02 0.0000000
# has-s         -0.0474128782 -0.06445150 -3.037426e-02 0.0000000
# is-s          -0.0387738877 -0.05188878 -2.565899e-02 0.0000000
# PL-GEN-s      -0.0385811893 -0.06170388 -1.545850e-02 0.0000211
# 3rdsg-plural  -0.0151149141 -0.03014508 -8.475015e-05 0.0476733
# GEN-plural    -0.0123355059 -0.02785690  3.185884e-03 0.2217609
# has-plural    -0.0263559131 -0.04506453 -7.647297e-03 0.0006941
# is-plural     -0.0177169225 -0.03293856 -2.495284e-03 0.0109241
# PL-GEN-plural -0.0175242241 -0.04190365  6.855203e-03 0.3384744
# GEN-3rdsg      0.0027794081 -0.01255425  1.811307e-02 0.9982993
# has-3rdsg     -0.0112409990 -0.02979416  7.312164e-03 0.5539156
# is-3rdsg      -0.0026020085 -0.01763217  1.242816e-02 0.9986883
# PL-GEN-3rdsg  -0.0024093101 -0.02666965  2.185103e-02 0.9999477
# has-GEN       -0.0140204071 -0.03297371  4.932893e-03 0.3034387
# is-GEN        -0.0053814166 -0.02090281  1.013997e-02 0.9480911
# PL-GEN-GEN    -0.0051887182 -0.02975642  1.937898e-02 0.9960127
# is-has         0.0086389905 -0.01006963  2.734761e-02 0.8198357
# PL-GEN-has     0.0088316889 -0.01786420  3.552758e-02 0.9585305
# PL-GEN-is      0.0001926984 -0.02418673  2.457213e-02 1.0000000
tuk <- TukeyHSD(aovmodel)
summary(tuk)
str(tuk)
tuk$morph
#                        diff         lwr           upr        p adj
# plural-s      -0.0210569652 -0.03417186 -7.942069e-03 5.174712e-05
# 3rdsg-s       -0.0361718792 -0.04906405 -2.327971e-02 0.000000e+00
# GEN-s         -0.0333924711 -0.04685411 -1.993083e-02 0.000000e+00
# has-s         -0.0474128782 -0.06445150 -3.037426e-02 0.000000e+00
# is-s          -0.0387738877 -0.05188878 -2.565899e-02 0.000000e+00
# PL-GEN-s      -0.0385811893 -0.06170388 -1.545850e-02 2.109079e-05
# 3rdsg-plural  -0.0151149141 -0.03014508 -8.475015e-05 4.767333e-02
# GEN-plural    -0.0123355059 -0.02785690  3.185884e-03 2.217609e-01
# has-plural    -0.0263559131 -0.04506453 -7.647297e-03 6.940932e-04
# is-plural     -0.0177169225 -0.03293856 -2.495284e-03 1.092407e-02
# PL-GEN-plural -0.0175242241 -0.04190365  6.855203e-03 3.384744e-01
# GEN-3rdsg      0.0027794081 -0.01255425  1.811307e-02 9.982993e-01
# has-3rdsg     -0.0112409990 -0.02979416  7.312164e-03 5.539156e-01
# is-3rdsg      -0.0026020085 -0.01763217  1.242816e-02 9.986883e-01
# PL-GEN-3rdsg  -0.0024093101 -0.02666965  2.185103e-02 9.999477e-01
# has-GEN       -0.0140204071 -0.03297371  4.932893e-03 3.034387e-01
# is-GEN        -0.0053814166 -0.02090281  1.013997e-02 9.480911e-01
# PL-GEN-GEN    -0.0051887182 -0.02975642  1.937898e-02 9.960127e-01
# is-has         0.0086389905 -0.01006963  2.734761e-02 8.198357e-01
# PL-GEN-has     0.0088316889 -0.01786420  3.552758e-02 9.585305e-01
# PL-GEN-is      0.0001926984 -0.02418673  2.457213e-02 1.000000e+00
bwplot(lenMorph ~ morph, d, ylab="duration of S (in sec)", cex.axis=1.5)
plot(lenMorph ~ morph, d, ylab="duration of S (in sec)", xlab="type of S")
dotplot(lenMorph ~ morph, d, ylab="duration of S (in sec)", cex.axis=1.5)
ScatterplotByFactor <- function (x, y, cex=0.25, lwd=2, linelength=0.75, linecolor="black", ylab=NULL, xlab=NULL, ylim=NULL, ...) {
if (is.null (ylim))
ylim <- range (x)
if (is.null (ylab))
ylab <- deparse (substitute (y))
if (is.null (xlab))
xlab <- deparse (substitute (x))
plot (NA, NA, type="n", xlim=c(0.5, nlevels (y) + 0.5), ylim=ylim, xaxt="n", ylab=ylab, xlab=xlab, ...)
points (jitter (as.numeric (y), 1.5), x, cex=cex, ...)
axis (1, 1:nlevels (y), levels (y))
means <- tapply (x, y, mean)
for (i in 1:nlevels (y))
lines (i + c (-0.5, +0.5) * linelength, c (means [i], means [i]), lwd=lwd, col=linecolor)
}
ScatterplotByFactor (d$lenMorph, d$morph, ylab="duration of S (sec)", xlab="type of S", las=1, lwd=3)
png ("scatterplot_by_typeOfS.png", width=850, height=550, pointsize=17)
ScatterplotByFactor (d$lenMorph, d$morph, ylab="duration of S (ms)", xlab="type of S", las=1, lwd=3)
dev.off()
pdf ("scatterplot_by_typeOfS.pdf", width=8.50, height=5.50, pointsize=10)
ScatterplotByFactor (d$lenMorph, d$morph, ylab="duration of S (ms)", xlab="type of S", las=1, lwd=3)
dev.off()
# this is our data set, which we will use for fitting the LMER models and the beta regressions:
# write.table(d, "data set prepared for analysis.csv", sep=",")
# we use a new method for comparing means in unbalanced designs
library(multcomp)
library(sandwich)
#create data with only two columns
names(d)
aovdat <- d[ , c(4,15)]
head(aovdat)
aovmodel <- aov(lenMorph ~ morph, aovdat)
aovmodel_glht <- glht(aovmodel, mcp(morph = "Tukey"), vcov = vcovHC)
summary(aovmodel_glht)
aovdat <- d[ , c(4,15)]
head(aovdat)
aovmodel <- aov(lenMorph ~ morph, aovdat)
aovmodel_glht <- glht(aovmodel, mcp(morph = "Tukey"), vcov = vcovHC)
summary(aovmodel_glht)
aovdat <- d[ , c(4,15)]
head(aovdat)
aovmodel <- aov(lenMorph ~ morph, aovdat)
aovmodel_glht <- glht(aovmodel, mcp(morph = "Tukey"), vcov = vcovHC)
summary(aovmodel_glht)
aovdat <- d[ , c(4,15)]
head(aovdat)
aovmodel <- aov(lenMorph ~ morph, aovdat)
aovmodel_glht <- glht(aovmodel, mcp(morph = "Tukey"), vcov = vcovHC)
summary(aovmodel_glht)
aovdat <- d[ , c(4,15)]
head(aovdat)
aovmodel <- aov(lenMorph ~ morph, aovdat)
aovmodel_glht <- glht(aovmodel, mcp(morph = "Tukey"), vcov = vcovHC)
summary(aovmodel_glht)
