#running R version 3.5.1 #Code written by Dani Crain, ddiancrain@gmail.com library(tidyverse) library(nlme) ICprogpreg=read.csv("IceSealClaws_prog_updatepregstats_W0only_Season_noNAs.csv", check.names=F) model_Z = aov(ProgZscores ~ Stage, data=ICprogpreg) anova(model_Z) #Analysis of Variance Table # #Response: ProgZscores #Df Sum Sq Mean Sq F value Pr(>F) #Stage 3 53.864 17.9548 6.5349 0.001928 ** # Residuals 26 71.436 2.7475 #--- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 ggplot(ICprogpreg, aes(x=Stage, y=ProgZscores))+ geom_boxplot() TukeyHSD(model_Z) #Tukey multiple comparisons of means #95% family-wise confidence level #Fit: aov(formula = ProgZscores ~ Stage, data = ICprogpreg) #$`Stage` #diff lwr upr p adj #late-early -0.6597635 -4.113637 2.7941103 0.9525150 #mid-early -2.0068130 -5.561551 1.5479250 0.4243068 #nonpreg-early -4.0512469 -7.764050 -0.3384442 0.0286626 #mid-late -1.3470495 -3.318863 0.6247644 0.2634489 #nonpreg-late -3.3914834 -5.635763 -1.1472036 0.0017105 #nonpreg-mid -2.0444340 -4.441038 0.3521699 0.1147821 plot(model_Z, 1) #install.packages("car") library(car) leveneTest(ProgZscores ~ Stage, data = ICprogpreg) #Levene's Test for Homogeneity of Variance (center = median) # Df F value Pr(>F) #group 3 0.8722 0.4668 # 29 png("PregStage_iceseals_noearly.png", width = 6, height = 4, units = 'in', res = 300) ggplot(ICprogpreg, aes(x=Stage, y=ProgZscores))+ geom_boxplot()+ geom_jitter(shape=21, position=position_jitter(0.2),size=2, aes(fill = Species))+ scale_x_discrete(limits=c("nonpreg", "mid","late"), labels=c("Non-pregnant\n(Fall & Winter)","Mid-pregnancy\n(Fall & Winter)","Late pregnancy\n(Spring & Summer)"))+ theme_classic()+ scale_fill_manual(values=c("gray","black"), labels=c("Bearded seal","Ringed seal"))+ xlab("Reproductive Stage")+ ylab("Progesterone Z-scores") dev.off()