TAWpupalheatstressanalysisCRM2026
Cailyn McKay
2026-03-27
Statistical Analyses for true armyworm (Mythimna unipuncta) pupal heat stress
project completed in 2022.
Armyworms were either heat stressed (30ºC) for: -14 hours (Treatment = 1) beginning on: prepupa (PP) day 1, pupa (P) day 1, P day 2, P day 3, or P day 4; correspond with Period levels 1-5 respectively
-48 hours (Treatment = 2) beginning on: PP day 1, P days 1-2, P days 3-4, P days 5-6, P days 7-8, or P days 9-10; correspond with Periods 6-11).
Controls were maintained at regular rearing conditions throughout development (25ºC, 60%RH, 16L:8D).
Cages consisted of eight males and five females and were composed of one of the following (variable name = “Cage_Setup” or “PopType”): control males and control females (Cage_Setup = 0), control males and treated females (Cage_Setup = 1), treated males and control females (Cage_Setup = 2), or treated males with treated females (Cage_Setup = 3; period treated always matched between the sexes i.e. Period 6 and Period 6)
There were 3 replicates of each. For each replicate I counted the number of total spermatophores (matings; Total_Spm) and calculated the number of spermatophores per emerged female (as some never emerged; Spm_perFemale). I also counted the number of eggs laid (Total_Eggs) and the number of those which were fertile (Fertile_Eggs) in each cage.
Importing and properly labelling factor data.
require(multcomp)
require(tidyverse)
require(Rmisc)
require(dplyr)
library(car)
library(readxl)
#import dataset. My code for my file...change as needed
TAW_pupalheatexp_14h_48h <- read_excel("TAW_pupalheatexp_14h_48h.xlsx",
sheet = "14h_48h_alldataforR")
#View(TAW_pupalheatexp_14h_48h)
str(TAW_pupalheatexp_14h_48h)Here I assign relevant variables as factors
TAW_pupalheatexp_14h_48h$Period <- factor(TAW_pupalheatexp_14h_48h$Period)
TAW_pupalheatexp_14h_48h$Treatment <- factor(TAW_pupalheatexp_14h_48h$Treatment)
#TAW_pupalheatexp_14h_48h$Cage_Setup <- factor(TAW_pupalheatexp_14h_48h$Cage_Setup)
#TAW_pupalheatexp_14h_48h$Female <- factor(TAW_pupalheatexp_14h_48h$Female)
#TAW_pupalheatexp_14h_48h$Male <- factor(TAW_pupalheatexp_14h_48h$Male)
TAW_pupalheatexp_14h_48h$Replicate <- factor(TAW_pupalheatexp_14h_48h$Replicate)
# to check call: str(TAW_pupalheatexp_14h_48h)Not shown: visualizing the dataset with boxplots for relevant dependent variables Looking at the boxplots, most of my results are seen in the 48 hour treatments so I’m dropping the 14h data.
TAW_48h <- read_excel("TAW_pupalheatexp_14h_48h.xlsx",
sheet = "48h_only")
#View(TAW_48h)Adding column for male treated and female treated (binary yes (1) or no (0)) first for full dataset then for 48h dataset.
I also define new variables here (PopType and AgeHS) to make my graphs simpler because they have character names.
Finally, ensuring my factors are labelled as factors. Also assigning new variables to english description rather than coded values.
Testing interactions between age at heat stress and the population type (or which sex(es) was treated)
Looking at Period to best identify the treatment and development stage which was most affected by heat stress. First identifying whether there are significant interactions due to sex which involves dropping the controls (unbalanced design) then continuing with single factor/main effects if needed.
Total number of spermatophores per population. Does not account for the (variable) number of eclosed females.
First removing controls from the dataset.
TAW_noCon <- subset(TAW_48h, Period !="0")Total spermatophores with sex interactions
SPMint <- lm(Spm_Total~Period *Cage_Setup, data=TAW_noCon)
summary(SPMint)
##
## Call:
## lm(formula = Spm_Total ~ Period * Cage_Setup, data = TAW_noCon)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.000 -1.667 0.000 1.917 4.667
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.3333 1.8509 2.341 0.0249 *
## Period7 2.3333 2.6176 0.891 0.3786
## Period8 -0.6667 2.6176 -0.255 0.8004
## Period9 3.0000 2.6176 1.146 0.2593
## Period10 2.3333 2.6176 0.891 0.3786
## Period11 3.3333 2.6176 1.273 0.2110
## Cage_Setup2 -1.3333 2.6176 -0.509 0.6136
## Cage_Setup3 0.6667 2.6176 0.255 0.8004
## Period7:Cage_Setup2 0.3333 3.7019 0.090 0.9288
## Period8:Cage_Setup2 5.0000 3.7019 1.351 0.1852
## Period9:Cage_Setup2 2.3333 3.7019 0.630 0.5325
## Period10:Cage_Setup2 3.0000 3.7019 0.810 0.4230
## Period11:Cage_Setup2 1.6667 3.7019 0.450 0.6552
## Period7:Cage_Setup3 0.3333 3.7019 0.090 0.9288
## Period8:Cage_Setup3 3.3333 3.7019 0.900 0.3739
## Period9:Cage_Setup3 3.6667 3.7019 0.990 0.3285
## Period10:Cage_Setup3 2.3333 3.7019 0.630 0.5325
## Period11:Cage_Setup3 -1.6667 3.7019 -0.450 0.6552
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.206 on 36 degrees of freedom
## Multiple R-squared: 0.3791, Adjusted R-squared: 0.08592
## F-statistic: 1.293 on 17 and 36 DF, p-value: 0.2514
anova(SPMint)
## Analysis of Variance Table
##
## Response: Spm_Total
## Df Sum Sq Mean Sq F value Pr(>F)
## Period 5 136.81 27.3630 2.6623 0.03794 *
## Cage_Setup 2 36.93 18.4630 1.7964 0.18045
## Period:Cage_Setup 10 52.19 5.2185 0.5077 0.87325
## Residuals 36 370.00 10.2778
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(SPMint)
Spermatophores per female with interactions. None significant
SPMperint <- lm(Spm_perFemale~Period *Cage_Setup, data=TAW_noCon)
summary(SPMperint)##
## Call:
## lm(formula = Spm_perFemale ~ Period * Cage_Setup, data = TAW_noCon)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.05000 -0.40000 0.02167 0.38333 1.31000
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.93333 0.37867 2.465 0.0186 *
## Period7 0.60000 0.53552 1.120 0.2700
## Period8 -0.13333 0.53552 -0.249 0.8048
## Period9 0.53333 0.53552 0.996 0.3259
## Period10 0.51667 0.53552 0.965 0.3411
## Period11 0.82333 0.53552 1.537 0.1329
## Cage_Setup2 0.25667 0.53552 0.479 0.6346
## Cage_Setup3 0.23333 0.53552 0.436 0.6656
## Period7:Cage_Setup2 -0.56667 0.75734 -0.748 0.4592
## Period8:Cage_Setup2 0.46000 0.75734 0.607 0.5474
## Period9:Cage_Setup2 -0.05667 0.75734 -0.075 0.9408
## Period10:Cage_Setup2 -0.04000 0.75734 -0.053 0.9582
## Period11:Cage_Setup2 -0.11333 0.75734 -0.150 0.8819
## Period7:Cage_Setup3 -0.16667 0.75734 -0.220 0.8271
## Period8:Cage_Setup3 0.50000 0.75734 0.660 0.5133
## Period9:Cage_Setup3 1.12333 0.75734 1.483 0.1467
## Period10:Cage_Setup3 0.25000 0.75734 0.330 0.7432
## Period11:Cage_Setup3 -0.59000 0.75734 -0.779 0.4410
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6559 on 36 degrees of freedom
## Multiple R-squared: 0.3917, Adjusted R-squared: 0.1044
## F-statistic: 1.364 on 17 and 36 DF, p-value: 0.2118anova(SPMperint)## Analysis of Variance Table
##
## Response: Spm_perFemale
## Df Sum Sq Mean Sq F value Pr(>F)
## Period 5 4.5891 0.91783 2.1336 0.08362 .
## Cage_Setup 2 1.5838 0.79191 1.8409 0.17330
## Period:Cage_Setup 10 3.7990 0.37990 0.8831 0.55720
## Residuals 36 15.4864 0.43018
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(SPMperint)
Total Eggs, also non-significant interactions
TEint <- lm(Total_Eggs~Period *Cage_Setup, data=TAW_noCon)
summary(TEint)##
## Call:
## lm(formula = Total_Eggs ~ Period * Cage_Setup, data = TAW_noCon)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2750.33 -670.92 34.33 742.08 2260.33
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1954.7 819.1 2.386 0.0224 *
## Period7 632.0 1158.4 0.546 0.5887
## Period8 -311.0 1158.4 -0.268 0.7899
## Period9 892.0 1158.4 0.770 0.4463
## Period10 1001.3 1158.4 0.864 0.3931
## Period11 656.7 1158.4 0.567 0.5743
## Cage_Setup2 -883.7 1158.4 -0.763 0.4505
## Cage_Setup3 -597.0 1158.4 -0.515 0.6094
## Period7:Cage_Setup2 474.0 1638.2 0.289 0.7740
## Period8:Cage_Setup2 1788.3 1638.2 1.092 0.2822
## Period9:Cage_Setup2 1875.0 1638.2 1.145 0.2599
## Period10:Cage_Setup2 1843.0 1638.2 1.125 0.2680
## Period11:Cage_Setup2 1234.3 1638.2 0.753 0.4561
## Period7:Cage_Setup3 -241.3 1638.2 -0.147 0.8837
## Period8:Cage_Setup3 215.3 1638.2 0.131 0.8962
## Period9:Cage_Setup3 1663.7 1638.2 1.016 0.3166
## Period10:Cage_Setup3 879.3 1638.2 0.537 0.5947
## Period11:Cage_Setup3 1467.0 1638.2 0.896 0.3765
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1419 on 36 degrees of freedom
## Multiple R-squared: 0.3731, Adjusted R-squared: 0.07705
## F-statistic: 1.26 on 17 and 36 DF, p-value: 0.2718anova(TEint)## Analysis of Variance Table
##
## Response: Total_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## Period 5 33494857 6698971 3.3283 0.01426 *
## Cage_Setup 2 1017002 508501 0.2526 0.77811
## Period:Cage_Setup 10 8609711 860971 0.4278 0.92335
## Residuals 36 72458748 2012743
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(TEint)
Fertile Eggs, interaction also non sig but model fit is also not
ideal.
FEint <- lm(Fertile_Eggs~Period *Cage_Setup, data=TAW_noCon)
summary(FEint)##
## Call:
## lm(formula = Fertile_Eggs ~ Period * Cage_Setup, data = TAW_noCon)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2642.67 -552.67 51.67 535.83 2130.33
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1533.3 750.2 2.044 0.0483 *
## Period7 552.0 1061.0 0.520 0.6061
## Period8 -328.3 1061.0 -0.309 0.7588
## Period9 150.7 1061.0 0.142 0.8879
## Period10 1063.3 1061.0 1.002 0.3229
## Period11 784.3 1061.0 0.739 0.4645
## Cage_Setup2 -632.7 1061.0 -0.596 0.5547
## Cage_Setup3 -340.3 1061.0 -0.321 0.7502
## Period7:Cage_Setup2 341.0 1500.4 0.227 0.8215
## Period8:Cage_Setup2 1405.0 1500.4 0.936 0.3553
## Period9:Cage_Setup2 2014.3 1500.4 1.342 0.1878
## Period10:Cage_Setup2 1289.3 1500.4 0.859 0.3959
## Period11:Cage_Setup2 1105.7 1500.4 0.737 0.4660
## Period7:Cage_Setup3 -832.0 1500.4 -0.555 0.5827
## Period8:Cage_Setup3 -181.7 1500.4 -0.121 0.9043
## Period9:Cage_Setup3 1809.3 1500.4 1.206 0.2357
## Period10:Cage_Setup3 803.3 1500.4 0.535 0.5957
## Period11:Cage_Setup3 923.3 1500.4 0.615 0.5422
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1299 on 36 degrees of freedom
## Multiple R-squared: 0.3864, Adjusted R-squared: 0.09666
## F-statistic: 1.334 on 17 and 36 DF, p-value: 0.2279anova(FEint)## Analysis of Variance Table
##
## Response: Fertile_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## Period 5 27565827 5513165 3.2651 0.01562 *
## Cage_Setup 2 1554636 777318 0.4604 0.63472
## Period:Cage_Setup 10 9159603 915960 0.5425 0.84817
## Residuals 36 60786026 1688501
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(FEint)
#### Here I am looking at whether the number of eggs covaries with the
number of spermatophores and testing the interactions of my other
factors as above
TEintANCOVA <- lm(Total_Eggs~Spm_Total*Period *Cage_Setup, data=TAW_noCon)
summary(TEintANCOVA)##
## Call:
## lm(formula = Total_Eggs ~ Spm_Total * Period * Cage_Setup, data = TAW_noCon)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1936.84 -276.93 -37.05 211.94 2260.35
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 261.97 1089.29 0.240 0.8125
## Spm_Total 390.62 191.57 2.039 0.0556 .
## Period7 2602.32 3987.10 0.653 0.5218
## Period8 577.23 1699.36 0.340 0.7378
## Period9 2584.58 2359.52 1.095 0.2870
## Period10 -143.51 1820.71 -0.079 0.9380
## Period11 1459.21 2526.49 0.578 0.5703
## Cage_Setup2 789.53 1980.54 0.399 0.6946
## Cage_Setup3 -4044.38 2551.00 -1.585 0.1294
## Spm_Total:Period7 -432.27 597.07 -0.724 0.4779
## Spm_Total:Period8 -171.22 355.30 -0.482 0.6354
## Spm_Total:Period9 -390.61 330.01 -1.184 0.2512
## Spm_Total:Period10 35.01 270.91 0.129 0.8985
## Spm_Total:Period11 -274.52 341.53 -0.804 0.4315
## Spm_Total:Cage_Setup2 -384.12 534.25 -0.719 0.4809
## Spm_Total:Cage_Setup3 637.39 468.89 1.359 0.1899
## Period7:Cage_Setup2 -3972.70 4538.43 -0.875 0.3923
## Period8:Cage_Setup2 -1596.18 2928.84 -0.545 0.5921
## Period9:Cage_Setup2 -627.03 4120.99 -0.152 0.8807
## Period10:Cage_Setup2 5410.55 4301.91 1.258 0.2237
## Period11:Cage_Setup2 -1241.85 3448.28 -0.360 0.7227
## Period7:Cage_Setup3 1547.19 5095.08 0.304 0.7647
## Period8:Cage_Setup3 2888.34 3598.72 0.803 0.4321
## Period9:Cage_Setup3 4243.79 5901.15 0.719 0.4808
## Period10:Cage_Setup3 2116.08 3660.94 0.578 0.5700
## Period11:Cage_Setup3 781.20 1721.42 0.454 0.6551
## Spm_Total:Period7:Cage_Setup2 866.21 806.78 1.074 0.2964
## Spm_Total:Period8:Cage_Setup2 507.78 648.58 0.783 0.4433
## Spm_Total:Period9:Cage_Setup2 483.58 689.53 0.701 0.4916
## Spm_Total:Period10:Cage_Setup2 -329.89 703.08 -0.469 0.6443
## Spm_Total:Period11:Cage_Setup2 479.66 633.17 0.758 0.4580
## Spm_Total:Period7:Cage_Setup3 -415.59 782.21 -0.531 0.6014
## Spm_Total:Period8:Cage_Setup3 -650.86 617.74 -1.054 0.3053
## Spm_Total:Period9:Cage_Setup3 -563.06 681.37 -0.826 0.4189
## Spm_Total:Period10:Cage_Setup3 -540.80 549.76 -0.984 0.3376
## Spm_Total:Period11:Cage_Setup3 NA NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1222 on 19 degrees of freedom
## Multiple R-squared: 0.7547, Adjusted R-squared: 0.3157
## F-statistic: 1.719 on 34 and 19 DF, p-value: 0.1063anova(TEintANCOVA)## Analysis of Variance Table
##
## Response: Total_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## Spm_Total 1 52797063 52797063 35.3783 1.003e-05 ***
## Period 5 8618415 1723683 1.1550 0.3667
## Cage_Setup 2 3322391 1661196 1.1131 0.3490
## Spm_Total:Period 5 2712595 542519 0.3635 0.8671
## Spm_Total:Cage_Setup 2 1236143 618072 0.4142 0.6667
## Period:Cage_Setup 10 8729335 872933 0.5849 0.8066
## Spm_Total:Period:Cage_Setup 9 9809565 1089952 0.7304 0.6770
## Residuals 19 28354811 1492358
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(TEintANCOVA, type = "II")## Note: model has aliased coefficients
## sums of squares computed by model comparison## Anova Table (Type II tests)
##
## Response: Total_Eggs
## Sum Sq Df F value Pr(>F)
## Spm_Total 30269096 1 20.2827 0.0002432 ***
## Period 7604322 5 1.0191 0.4340425
## Cage_Setup 3068641 2 1.0281 0.3767404
## Spm_Total:Period 3270494 5 0.4383 0.8162525
## Spm_Total:Cage_Setup 520572 2 0.1744 0.8412799
## Period:Cage_Setup 8729335 10 0.5849 0.8066280
## Spm_Total:Period:Cage_Setup 9809565 9 0.7304 0.6769794
## Residuals 28354811 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(TEintANCOVA)## Warning: not plotting observations with leverage one:
## 6, 21, 54
## Warning in sqrt(crit * p * (1 - hh)/hh): NaNs produced## Warning in sqrt(crit * p * (1 - hh)/hh): NaNs produced
Fertile Eggs, same as above…not pretty but not interacting. Can continue with main effects only now including controls.
FEintANCOVA <- lm(Fertile_Eggs~Spm_Total*Period *Cage_Setup, data=TAW_noCon)
summary(FEintANCOVA)##
## Call:
## lm(formula = Fertile_Eggs ~ Spm_Total * Period * Cage_Setup,
## data = TAW_noCon)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1778.42 -197.55 -5.76 214.02 1397.96
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -62.43 855.99 -0.073 0.9426
## Spm_Total 368.25 150.54 2.446 0.0243 *
## Period7 2195.86 3133.16 0.701 0.4919
## Period8 113.53 1335.40 0.085 0.9331
## Period9 923.56 1854.17 0.498 0.6241
## Period10 -471.99 1430.76 -0.330 0.7451
## Period11 1549.36 1985.38 0.780 0.4448
## Cage_Setup2 757.10 1556.36 0.486 0.6322
## Cage_Setup3 -3416.55 2004.64 -1.704 0.1046
## Spm_Total:Period7 -375.47 469.19 -0.800 0.4335
## Spm_Total:Period8 -53.55 279.20 -0.192 0.8499
## Spm_Total:Period9 -256.04 259.33 -0.987 0.3359
## Spm_Total:Period10 101.41 212.89 0.476 0.6393
## Spm_Total:Period11 -259.90 268.39 -0.968 0.3450
## Spm_Total:Cage_Setup2 -299.59 419.83 -0.714 0.4842
## Spm_Total:Cage_Setup3 566.14 368.46 1.536 0.1409
## Period7:Cage_Setup2 -3278.88 3566.41 -0.919 0.3694
## Period8:Cage_Setup2 -2027.75 2301.55 -0.881 0.3893
## Period9:Cage_Setup2 -383.70 3238.38 -0.118 0.9069
## Period10:Cage_Setup2 2937.71 3380.55 0.869 0.3957
## Period11:Cage_Setup2 -918.32 2709.75 -0.339 0.7384
## Period7:Cage_Setup3 921.11 4003.84 0.230 0.8205
## Period8:Cage_Setup3 2830.19 2827.96 1.001 0.3295
## Period9:Cage_Setup3 5300.54 4637.27 1.143 0.2672
## Period10:Cage_Setup3 2347.07 2876.86 0.816 0.4247
## Period11:Cage_Setup3 333.62 1352.73 0.247 0.8078
## Spm_Total:Period7:Cage_Setup2 691.86 633.99 1.091 0.2888
## Spm_Total:Period8:Cage_Setup2 420.83 509.67 0.826 0.4192
## Spm_Total:Period9:Cage_Setup2 407.11 541.85 0.751 0.4617
## Spm_Total:Period10:Cage_Setup2 -158.92 552.49 -0.288 0.7767
## Spm_Total:Period11:Cage_Setup2 374.35 497.56 0.752 0.4610
## Spm_Total:Period7:Cage_Setup3 -392.62 614.68 -0.639 0.5306
## Spm_Total:Period8:Cage_Setup3 -721.94 485.44 -1.487 0.1534
## Spm_Total:Period9:Cage_Setup3 -643.39 535.44 -1.202 0.2443
## Spm_Total:Period10:Cage_Setup3 -553.37 432.02 -1.281 0.2156
## Spm_Total:Period11:Cage_Setup3 NA NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 960 on 19 degrees of freedom
## Multiple R-squared: 0.8233, Adjusted R-squared: 0.507
## F-statistic: 2.603 on 34 and 19 DF, p-value: 0.01513anova(FEintANCOVA)## Analysis of Variance Table
##
## Response: Fertile_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## Spm_Total 1 49786930 49786930 54.0245 5.761e-07 ***
## Period 5 9129450 1825890 1.9813 0.1278
## Cage_Setup 2 4146976 2073488 2.2500 0.1327
## Spm_Total:Period 5 2679096 535819 0.5814 0.7139
## Spm_Total:Cage_Setup 2 878828 439414 0.4768 0.6280
## Period:Cage_Setup 10 8163742 816374 0.8859 0.5621
## Spm_Total:Period:Cage_Setup 9 6771405 752378 0.8164 0.6082
## Residuals 19 17509664 921561
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(FEintANCOVA)## Warning: not plotting observations with leverage one:
## 6, 21, 54
## Warning in sqrt(crit * p * (1 - hh)/hh): NaNs produced## Warning in sqrt(crit * p * (1 - hh)/hh): NaNs produced
Now my models will focus on main effects only. In all cases I use type II ANOVAs due to the unbalanced design of the study. For simplicity, I will be dropping the 14h treatments.
Matings (via Spermatophore counts)
Starting with the total number of spermatophores per population
#spm total
SPMtot <- lm(Spm_Total~AgeHS + Male + Female, data=TAW_48h)
summary(SPMtot) ##
## Call:
## lm(formula = Spm_Total ~ AgeHS + Male + Female, data = TAW_48h)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.3148 -2.0093 -0.2037 1.9074 5.7963
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 12.5000 1.4763 8.467 3.72e-11 ***
## AgeHSPrepupa -10.5741 2.1071 -5.018 7.26e-06 ***
## AgeHSDays 1-2 -8.0185 2.1071 -3.806 0.000394 ***
## AgeHSDays 3-4 -8.4630 2.1071 -4.016 0.000203 ***
## AgeHSDays 5-6 -5.5741 2.1071 -2.645 0.010935 *
## AgeHSDays 7-8 -6.4630 2.1071 -3.067 0.003513 **
## AgeHSDays 9-10 -7.2407 2.1071 -3.436 0.001210 **
## Male1 2.0000 0.9842 2.032 0.047585 *
## Female1 1.2778 0.9842 1.298 0.200271
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.953 on 49 degrees of freedom
## Multiple R-squared: 0.4026, Adjusted R-squared: 0.3051
## F-statistic: 4.128 on 8 and 49 DF, p-value: 0.0007964anova(SPMtot)## Analysis of Variance Table
##
## Response: Spm_Total
## Df Sum Sq Mean Sq F value Pr(>F)
## AgeHS 6 250.99 41.832 4.7983 0.0006387 ***
## Male 1 22.23 22.231 2.5500 0.1167208
## Female 1 14.69 14.694 1.6855 0.2002706
## Residuals 49 427.19 8.718
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(SPMtot, type = "II")## Anova Table (Type II tests)
##
## Response: Spm_Total
## Sum Sq Df F value Pr(>F)
## AgeHS 279.93 6 5.3516 0.0002589 ***
## Male 36.00 1 4.1294 0.0475850 *
## Female 14.69 1 1.6855 0.2002706
## Residuals 427.19 49
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(SPMtot)
SPMsumPeriod <- summarySE(data = TAW_48h, measurevar = 'Spm_Total', groupvars = 'AgeHS' )
SPMsumPeriod## AgeHS N Spm_Total sd se ci
## 1 Control 4 12.500000 1.290994 0.6454972 2.054260
## 2 Prepupa 9 4.111111 2.571208 0.8570694 1.976406
## 3 Days 1-2 9 6.666667 2.828427 0.9428090 2.174122
## 4 Days 3-4 9 6.222222 3.527668 1.1758895 2.711606
## 5 Days 5-6 9 9.111111 3.018462 1.0061539 2.320195
## 6 Days 7-8 9 8.222222 3.456074 1.1520245 2.656573
## 7 Days 9-10 9 7.444444 3.045944 1.0153148 2.341320SPMsumSex <- summarySE(data = TAW_48h, measurevar = 'Spm_Total', groupvars = 'PopType' )
SPMsumSex## PopType N Spm_Total sd se ci
## 1 Neither Sex 4 12.500000 1.290994 0.6454972 2.054260
## 2 Males Only 18 6.777778 3.439486 0.8106947 1.710416
## 3 Both Sexes 18 8.055556 3.207935 0.7561176 1.595269
## 4 Females Only 18 6.055556 3.280463 0.7732126 1.631336Post hoc for period and male status
Now for Dunnett contrasts. Note, throughout this work I am only interested in the comparisons to the control, not between treatments.
library(emmeans)## Welcome to emmeans.
## Caution: You lose important information if you filter this package's results.
## See '? untidy'SPM.means <- emmeans(SPMtot, ~AgeHS)
SPM.Malemeans <- emmeans(SPMtot, ~Male)
# Dunnett contrasts
contrast(SPM.means, adjust="bonferroni", method="dunnett")## contrast estimate SE df t.ratio p.value
## Prepupa - Control -10.57 2.11 49 -5.018 <.0001
## (Days 1-2) - Control -8.02 2.11 49 -3.806 0.0024
## (Days 3-4) - Control -8.46 2.11 49 -4.016 0.0012
## (Days 5-6) - Control -5.57 2.11 49 -2.645 0.0656
## (Days 7-8) - Control -6.46 2.11 49 -3.067 0.0211
## (Days 9-10) - Control -7.24 2.11 49 -3.436 0.0073
##
## Results are averaged over the levels of: Male, Female
## P value adjustment: bonferroni method for 6 testscontrast(SPM.Malemeans, adjust="bonferroni", method="dunnett")## contrast estimate SE df t.ratio p.value
## Male1 - Male0 2 0.984 49 2.032 0.0476
##
## Results are averaged over the levels of: AgeHS, FemaleTreatments different from the controls (all by decreases):
Spermatophores per female.
Effectively measuring the average number of times a female mated in her lifetime. Each replicate has one value which is the total # spermatophores (as above) divided by the number of emerged females. Alternatively this corrects for the fact that in some cages, not all 5 females emerged so the above number is an unfair assessment of actual mating “rates”. As above, the model looks at the type of treatment through the period it occurred and what type of males she had available.
#spm per female
SPMperF <- lm(Spm_perFemale~AgeHS+Male+Female, data=TAW_48h)
summary(SPMperF) ##
## Call:
## lm(formula = Spm_perFemale ~ AgeHS + Male + Female, data = TAW_48h)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.07556 -0.47167 -0.02639 0.33319 1.47278
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.5000 0.3153 7.929 2.45e-10 ***
## AgeHSPrepupa -1.8267 0.4500 -4.059 0.000177 ***
## AgeHSDays 1-2 -1.4711 0.4500 -3.269 0.001977 **
## AgeHSDays 3-4 -1.6400 0.4500 -3.644 0.000647 ***
## AgeHSDays 5-6 -0.9378 0.4500 -2.084 0.042407 *
## AgeHSDays 7-8 -1.2400 0.4500 -2.755 0.008206 **
## AgeHSDays 9-10 -1.2378 0.4500 -2.751 0.008313 **
## Male1 0.4194 0.2102 1.995 0.051570 .
## Female1 0.2156 0.2102 1.025 0.310176
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6306 on 49 degrees of freedom
## Multiple R-squared: 0.3316, Adjusted R-squared: 0.2225
## F-statistic: 3.039 on 8 and 49 DF, p-value: 0.007437anova(SPMperF)## Analysis of Variance Table
##
## Response: Spm_perFemale
## Df Sum Sq Mean Sq F value Pr(>F)
## AgeHS 6 8.0852 1.34753 3.3886 0.007106 **
## Male 1 1.1656 1.16563 2.9312 0.093205 .
## Female 1 0.4182 0.41818 1.0516 0.310176
## Residuals 49 19.4854 0.39766
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(SPMperF, type = "II")## Anova Table (Type II tests)
##
## Response: Spm_perFemale
## Sum Sq Df F value Pr(>F)
## AgeHS 9.2410 6 3.8731 0.003055 **
## Male 1.5834 1 3.9818 0.051570 .
## Female 0.4182 1 1.0516 0.310176
## Residuals 19.4854 49
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(SPMperF)
#summaries data
SPMpersumPeriod <- summarySE(data = TAW_48h, measurevar = 'Spm_perFemale', groupvars = 'AgeHS' )
SPMpersumPeriod## AgeHS N Spm_perFemale sd se ci
## 1 Control 4 2.500000 0.2581989 0.1290994 0.4108521
## 2 Prepupa 9 1.096667 0.7416030 0.2472010 0.5700465
## 3 Days 1-2 9 1.452222 0.4864868 0.1621623 0.3739469
## 4 Days 3-4 9 1.283333 0.6595453 0.2198484 0.5069714
## 5 Days 5-6 9 1.985556 0.8621936 0.2873979 0.6627407
## 6 Days 7-8 9 1.683333 0.6855655 0.2285218 0.5269723
## 7 Days 9-10 9 1.685556 0.4166867 0.1388956 0.3202937SPMpersumSex <- summarySE(data = TAW_48h, measurevar = 'Spm_perFemale', groupvars = 'PopType' )
SPMpersumSex## PopType N Spm_perFemale sd se ci
## 1 Neither Sex 4 2.500000 0.2581989 0.1290994 0.4108521
## 2 Males Only 18 1.527222 0.6442945 0.1518617 0.3204001
## 3 Both Sexes 18 1.742778 0.7484478 0.1764108 0.3721944
## 4 Females Only 18 1.323333 0.6550528 0.1543974 0.3257501Period: F(6,49) = 3.8731, p < 0.003** Male: F(1,49) = 3.9818, p =
0.0516 . Female: F(1,49) = 1.0516, p = 0.3102
Residuals - 19.4854
Once again, the period is a significant predictor of the number of matings. Moving to Dunnett contrasts:
SPMper.means <- emmeans(SPMperF, ~AgeHS)
SPMper.Malemeans <- emmeans(SPMperF, ~Male)
# Dunnett contrasts
contrast(SPMper.means, adjust="bonferroni", method="dunnett")## contrast estimate SE df t.ratio p.value
## Prepupa - Control -1.827 0.45 49 -4.059 0.0011
## (Days 1-2) - Control -1.471 0.45 49 -3.269 0.0119
## (Days 3-4) - Control -1.640 0.45 49 -3.644 0.0039
## (Days 5-6) - Control -0.938 0.45 49 -2.084 0.2544
## (Days 7-8) - Control -1.240 0.45 49 -2.755 0.0492
## (Days 9-10) - Control -1.238 0.45 49 -2.751 0.0499
##
## Results are averaged over the levels of: Male, Female
## P value adjustment: bonferroni method for 6 testscontrast(SPMper.Malemeans, adjust="bonferroni", method="dunnett")## contrast estimate SE df t.ratio p.value
## Male1 - Male0 0.419 0.21 49 1.995 0.0516
##
## Results are averaged over the levels of: AgeHS, FemaleSignificantly different from the controls (all by decreases): Period 6: 48h beginning on first day as prepupa Period 7: 48h on days 1-2 Period 8: 48h on days 3-4 pupa Period 10: 48h on days 7-8 Period 11: 48h on days 9-10
Treated males once again had a positive effect (+ 0.419) on the number of matings which occurred when controlled for treatment.
Now models for the total number of eggs laid in the cage (population) over its lifespan.
ANOVA_TE <- lm(Total_Eggs~AgeHS + Male + Female ,data=TAW_48h)
summary(ANOVA_TE)##
## Call:
## lm(formula = Total_Eggs ~ AgeHS + Male + Female, data = TAW_48h)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2820.30 -799.71 61.15 682.43 2697.37
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5128.5 667.7 7.681 5.9e-10 ***
## AgeHSPrepupa -3544.2 953.0 -3.719 0.000515 ***
## AgeHSDays 1-2 -2834.7 953.0 -2.975 0.004544 **
## AgeHSDays 3-4 -3187.3 953.0 -3.345 0.001587 **
## AgeHSDays 5-6 -1472.7 953.0 -1.545 0.128704
## AgeHSDays 7-8 -1635.4 953.0 -1.716 0.092455 .
## AgeHSDays 9-10 -1987.1 953.0 -2.085 0.042287 *
## Male1 67.0 445.1 0.151 0.880975
## Female1 -251.8 445.1 -0.566 0.574231
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1335 on 49 degrees of freedom
## Multiple R-squared: 0.4032, Adjusted R-squared: 0.3058
## F-statistic: 4.139 on 8 and 49 DF, p-value: 0.0007796anova(ANOVA_TE)## Analysis of Variance Table
##
## Response: Total_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## AgeHS 6 58030066 9671678 5.4235 0.0002307 ***
## Male 1 446473 446473 0.2504 0.6190587
## Female 1 570528 570528 0.3199 0.5742312
## Residuals 49 87381468 1783295
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(ANOVA_TE, type = "II")## Anova Table (Type II tests)
##
## Response: Total_Eggs
## Sum Sq Df F value Pr(>F)
## AgeHS 47825170 6 4.4697 0.001105 **
## Male 40401 1 0.0227 0.880975
## Female 570528 1 0.3199 0.574231
## Residuals 87381468 49
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(ANOVA_TE)
TEsumPeriod <- summarySE(data = TAW_48h, measurevar = 'Total_Eggs', groupvars = 'AgeHS' )
TEsumPeriod## AgeHS N Total_Eggs sd se ci
## 1 Control 4 5128.500 1450.633 725.3165 2308.2809
## 2 Prepupa 9 1461.111 1082.670 360.8900 832.2138
## 3 Days 1-2 9 2170.667 1060.018 353.3392 814.8017
## 4 Days 3-4 9 1818.000 1203.939 401.3129 925.4292
## 5 Days 5-6 9 3532.667 1521.323 507.1076 1169.3923
## 6 Days 7-8 9 3369.889 1634.647 544.8825 1256.5012
## 7 Days 9-10 9 3018.222 1236.490 412.1633 950.4502Adj R-squared = 0.31 Period: F(6,49) = 4.4697, p = 0.001** Female:
F(1,49) = 0.3199, p = 0.5742
Male: F(1,49) = 0.0227, p = 0.8810
For the total number of eggs laid (fecundity) the period and thus treatment were significant but regardless of which sex was treated. Dunnett contrasts:
TE.means <- emmeans(ANOVA_TE, ~AgeHS)
# Dunnett contrasts
contrast(TE.means, adjust="bonferroni", method="dunnett")## contrast estimate SE df t.ratio p.value
## Prepupa - Control -3544 953 49 -3.719 0.0031
## (Days 1-2) - Control -2835 953 49 -2.975 0.0273
## (Days 3-4) - Control -3187 953 49 -3.345 0.0095
## (Days 5-6) - Control -1473 953 49 -1.545 0.7722
## (Days 7-8) - Control -1635 953 49 -1.716 0.5547
## (Days 9-10) - Control -1987 953 49 -2.085 0.2537
##
## Results are averaged over the levels of: Male, Female
## P value adjustment: bonferroni method for 6 testsSignificant differences from controls: Period 6 Period 7 Period 8 Both we saw the have significant differences in mating as well…
Conclusion: Fecundity significantly affected by pupal heat stress, especially in early treatment periods of 48h duration.
Effect of pupal heat stress on number of fertile eggs:
ANOVA_FE <- lm(Fertile_Eggs~AgeHS + Male+Female,data=TAW_48h)
summary(ANOVA_FE)##
## Call:
## lm(formula = Fertile_Eggs ~ AgeHS + Male + Female, data = TAW_48h)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2551.13 -711.84 87.45 633.94 2439.04
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4929.50 628.22 7.847 3.28e-10 ***
## AgeHSPrepupa -3565.09 896.63 -3.976 0.000230 ***
## AgeHSDays 1-2 -3176.76 896.63 -3.543 0.000880 ***
## AgeHSDays 3-4 -3485.65 896.63 -3.887 0.000305 ***
## AgeHSDays 5-6 -2139.87 896.63 -2.387 0.020912 *
## AgeHSDays 7-8 -1804.20 896.63 -2.012 0.049712 *
## AgeHSDays 9-10 -2104.43 896.63 -2.347 0.023009 *
## Male1 80.06 418.82 0.191 0.849200
## Female1 -313.17 418.82 -0.748 0.458188
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1256 on 49 degrees of freedom
## Multiple R-squared: 0.4358, Adjusted R-squared: 0.3437
## F-statistic: 4.731 on 8 and 49 DF, p-value: 0.0002432anova(ANOVA_FE)## Analysis of Variance Table
##
## Response: Fertile_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## AgeHS 6 58200022 9700004 6.1445 7.423e-05 ***
## Male 1 671976 671976 0.4257 0.5172
## Female 1 882660 882660 0.5591 0.4582
## Residuals 49 77354186 1578657
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(ANOVA_FE, type = "II")## Anova Table (Type II tests)
##
## Response: Fertile_Eggs
## Sum Sq Df F value Pr(>F)
## AgeHS 45226372 6 4.7748 0.000664 ***
## Male 57680 1 0.0365 0.849200
## Female 882660 1 0.5591 0.458188
## Residuals 77354186 49
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(ANOVA_FE)
FEsumPeriod <- summarySE(data = TAW_48h, measurevar = 'Fertile_Eggs', groupvars = 'AgeHS' )
FEsumPeriod## AgeHS N Fertile_Eggs sd se ci
## 1 Control 4 4929.500 1571.470 785.7352 2500.5600
## 2 Prepupa 9 1209.000 1004.576 334.8586 772.1854
## 3 Days 1-2 9 1597.333 1012.680 337.5600 778.4146
## 4 Days 3-4 9 1288.444 1229.301 409.7669 944.9242
## 5 Days 5-6 9 2634.222 1245.382 415.1272 957.2851
## 6 Days 7-8 9 2969.889 1597.276 532.4253 1227.7748
## 7 Days 9-10 9 2669.667 1135.515 378.5050 872.8341Period: F(6,49) = 4.775, p = 0.007*** Male: F(1,49) = 0.037, p =
0.849
Female: F(1, 49) = 0.559, p = 0.458
Significant effect of treatment/period on the number of fertile eggs laid. Dunnett contrasts:
FE.means <- emmeans(ANOVA_FE, ~AgeHS)
# Dunnett contrasts
contrast(FE.means, adjust="bonferroni", method="dunnett")## contrast estimate SE df t.ratio p.value
## Prepupa - Control -3565 897 49 -3.976 0.0014
## (Days 1-2) - Control -3177 897 49 -3.543 0.0053
## (Days 3-4) - Control -3486 897 49 -3.887 0.0018
## (Days 5-6) - Control -2140 897 49 -2.387 0.1255
## (Days 7-8) - Control -1804 897 49 -2.012 0.2983
## (Days 9-10) - Control -2104 897 49 -2.347 0.1381
##
## Results are averaged over the levels of: Male, Female
## P value adjustment: bonferroni method for 6 testsSignificant periods: Period 6 Period 7 Period 8 Also very significant in similar treatment groups.
Now I am going to consider proportion fertile using a binomial glm with the logit link function. This calculates the probability of obtaining a value of 1 (fully fertile) for my dependent variable.
TAW_48h$PropFertile <- TAW_48h$Fertile_Eggs/TAW_48h$Total_Eggs
PropFEglm <- glm(PropFertile ~ AgeHS + Male +Female, data = TAW_48h, family = binomial, weight = Total_Eggs)
summary(PropFEglm)##
## Call:
## glm(formula = PropFertile ~ AgeHS + Male + Female, family = binomial,
## data = TAW_48h, weights = Total_Eggs)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.20969 0.03615 88.782 <2e-16 ***
## AgeHSPrepupa -1.38198 0.04721 -29.273 <2e-16 ***
## AgeHSDays 1-2 -1.95826 0.04393 -44.574 <2e-16 ***
## AgeHSDays 3-4 -2.14224 0.04406 -48.617 <2e-16 ***
## AgeHSDays 5-6 -1.91919 0.04351 -44.113 <2e-16 ***
## AgeHSDays 7-8 -0.99417 0.04493 -22.130 <2e-16 ***
## AgeHSDays 9-10 -0.94305 0.04585 -20.567 <2e-16 ***
## Male1 0.02061 0.01686 1.223 0.222
## Female1 -0.34892 0.01721 -20.271 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 28246 on 57 degrees of freedom
## Residual deviance: 19006 on 49 degrees of freedom
## AIC: 19442
##
## Number of Fisher Scoring iterations: 5car::Anova(PropFEglm, type = "II")## Analysis of Deviance Table (Type II tests)
##
## Response: PropFertile
## LR Chisq Df Pr(>Chisq)
## AgeHS 6814.3 6 <2e-16 ***
## Male 1.5 1 0.2215
## Female 412.5 1 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1pFEsumPeriod <- summarySE(data = TAW_48h, measurevar = 'PropFertile', groupvars = 'AgeHS' )
pFEsumPeriod## AgeHS N PropFertile sd se ci
## 1 Control 4 0.9514276 0.04940313 0.02470156 0.07861140
## 2 Prepupa 9 0.7678551 0.22032778 0.07344259 0.16935893
## 3 Days 1-2 9 0.6247654 0.28003373 0.09334458 0.21525298
## 4 Days 3-4 9 0.5403305 0.35060500 0.11686833 0.26949886
## 5 Days 5-6 9 0.7239127 0.17228735 0.05742912 0.13243178
## 6 Days 7-8 9 0.8307650 0.21960172 0.07320057 0.16880083
## 7 Days 9-10 9 0.8804201 0.08392655 0.02797552 0.06451165pFEsumSex <- summarySE(data = TAW_48h, measurevar = 'PropFertile', groupvars = 'PopType' )
pFEsumSex## PopType N PropFertile sd se ci
## 1 Neither Sex 4 0.9514276 0.04940313 0.02470156 0.0786114
## 2 Males Only 18 0.7636436 0.27413910 0.06461521 0.1363262
## 3 Both Sexes 18 0.7216121 0.23860138 0.05623888 0.1186537
## 4 Females Only 18 0.6987688 0.25733825 0.06065521 0.1279713PropFE.means.Period <- emmeans(PropFEglm, ~ AgeHS|Male+Female)
PropFEmeans.Female <- emmeans(PropFEglm, ~Female|AgeHS)
logit.dunnett <- contrast(emmeans(PropFEglm, ~ AgeHS), adjust ="bonferroni", method = "dunnett")
summary(logit.dunnett)## contrast estimate SE df z.ratio p.value
## Prepupa - Control -1.382 0.0472 Inf -29.273 <.0001
## (Days 1-2) - Control -1.958 0.0439 Inf -44.574 <.0001
## (Days 3-4) - Control -2.142 0.0441 Inf -48.617 <.0001
## (Days 5-6) - Control -1.919 0.0435 Inf -44.113 <.0001
## (Days 7-8) - Control -0.994 0.0449 Inf -22.130 <.0001
## (Days 9-10) - Control -0.943 0.0459 Inf -20.567 <.0001
##
## Results are averaged over the levels of: Male, Female
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: bonferroni method for 6 testslogit.dunnett.fem <- contrast(emmeans(PropFEglm, ~ Female, adjust ="bonferroni", method = "dunnett"))
summary(logit.dunnett.fem)## contrast estimate SE df z.ratio p.value
## Female0 effect 0.174 0.00861 Inf 20.271 <.0001
## Female1 effect -0.174 0.00861 Inf -20.271 <.0001
##
## Results are averaged over the levels of: AgeHS, Male
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: fdr method for 2 tests#backtransformed summary
summary(logit.dunnett, type = "response")## contrast odds.ratio SE df null z.ratio p.value
## Prepupa / Control 0.251 0.01190 Inf 1 -29.273 <.0001
## (Days 1-2) / Control 0.141 0.00620 Inf 1 -44.574 <.0001
## (Days 3-4) / Control 0.117 0.00517 Inf 1 -48.617 <.0001
## (Days 5-6) / Control 0.147 0.00638 Inf 1 -44.113 <.0001
## (Days 7-8) / Control 0.370 0.01660 Inf 1 -22.130 <.0001
## (Days 9-10) / Control 0.389 0.01790 Inf 1 -20.567 <.0001
##
## Results are averaged over the levels of: Male, Female
## P value adjustment: bonferroni method for 6 tests
## Tests are performed on the log odds ratio scaleAll treatment groups had a significantly reduced proportion of fertile
Known that multiple matings increases fecundity so accounting for the statistically significant differences in the incidence of mating following heat stress. Using ANCOVA, now including the number of spermatophores (successful matings), as the covariate. This is the ANCOVA considering the number of matings. Interestingly, we see the same trend if we include the 14h dataset.
SPMfactTE <- lm(Total_Eggs~Spm_Total + AgeHS + Male + Female,data=TAW_48h)
summary(SPMfactTE)##
## Call:
## lm(formula = Total_Eggs ~ Spm_Total + AgeHS + Male + Female,
## data = TAW_48h)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1920.37 -528.48 -41.66 504.36 2843.17
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1734.93 846.79 2.049 0.046 *
## Spm_Total 271.49 52.21 5.200 4.07e-06 ***
## AgeHSPrepupa -673.50 947.54 -0.711 0.481
## AgeHSDays 1-2 -657.74 876.53 -0.750 0.457
## AgeHSDays 3-4 -889.75 887.85 -1.002 0.321
## AgeHSDays 5-6 40.63 823.24 0.049 0.961
## AgeHSDays 7-8 119.17 840.77 0.142 0.888
## AgeHSDays 9-10 -21.34 857.88 -0.025 0.980
## Male1 -475.97 374.56 -1.271 0.210
## Female1 -598.68 365.84 -1.636 0.108
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1079 on 48 degrees of freedom
## Multiple R-squared: 0.6183, Adjusted R-squared: 0.5467
## F-statistic: 8.638 on 9 and 48 DF, p-value: 1.478e-07anova(SPMfactTE)## Analysis of Variance Table
##
## Response: Total_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## Spm_Total 1 76124065 76124065 65.3705 1.638e-10 ***
## AgeHS 6 11013922 1835654 1.5763 0.1746
## Male 1 276000 276000 0.2370 0.6286
## Female 1 3118443 3118443 2.6779 0.1083
## Residuals 48 55896106 1164502
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(SPMfactTE, type = "II")## Anova Table (Type II tests)
##
## Response: Total_Eggs
## Sum Sq Df F value Pr(>F)
## Spm_Total 31485362 1 27.0376 4.069e-06 ***
## AgeHS 7449198 6 1.0661 0.3958
## Male 1880460 1 1.6148 0.2099
## Female 3118443 1 2.6779 0.1083
## Residuals 55896106 48
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(SPMfactTE)
SPMfactFE <- lm(Fertile_Eggs~Spm_Total+AgeHS+Male +Female,data=TAW_48h)
summary(SPMfactFE)##
## Call:
## lm(formula = Fertile_Eggs ~ Spm_Total + AgeHS + Male + Female,
## data = TAW_48h)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1721.78 -520.28 -16.19 574.98 2593.60
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1331.80 733.74 1.815 0.0758 .
## Spm_Total 287.82 45.24 6.362 7e-08 ***
## AgeHSPrepupa -521.70 821.03 -0.635 0.5282
## AgeHSDays 1-2 -868.90 759.50 -1.144 0.2583
## AgeHSDays 3-4 -1049.87 769.31 -1.365 0.1787
## AgeHSDays 5-6 -535.56 713.33 -0.751 0.4564
## AgeHSDays 7-8 55.94 728.52 0.077 0.9391
## AgeHSDays 9-10 -20.42 743.34 -0.027 0.9782
## Male1 -495.58 324.55 -1.527 0.1333
## Female1 -680.93 317.00 -2.148 0.0368 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 935 on 48 degrees of freedom
## Multiple R-squared: 0.6939, Adjusted R-squared: 0.6365
## F-statistic: 12.09 on 9 and 48 DF, p-value: 1.077e-09anova(SPMfactFE)## Analysis of Variance Table
##
## Response: Fertile_Eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## Spm_Total 1 76984594 76984594 88.0516 1.949e-12 ***
## AgeHS 6 13926990 2321165 2.6548 0.02639 *
## Male 1 196038 196038 0.2242 0.63799
## Female 1 4034240 4034240 4.6142 0.03678 *
## Residuals 48 41966982 874312
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(SPMfactFE, type = "II")## Anova Table (Type II tests)
##
## Response: Fertile_Eggs
## Sum Sq Df F value Pr(>F)
## Spm_Total 35387204 1 40.4743 7.003e-08 ***
## AgeHS 8701068 6 1.6587 0.15181
## Male 2038570 1 2.3316 0.13333
## Female 4034240 1 4.6142 0.03678 *
## Residuals 41966982 48
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(SPMfactFE)###Below is my analysis for the pupal mortality rates. This is in a different dataset. Same thing, ensuring all factors in the dataset are treated as such.
TAW_pupalheatexp_pupalmort <- read_excel("TAW_pupalheatexp_14h_48h.xlsx",
sheet = "mort_all")
#View(TAW_pupalheatexp_pupalmort)
TAW_pupalheatexp_pupalmort$Period <- factor(TAW_pupalheatexp_pupalmort$Period)
TAW_pupalheatexp_pupalmort$Treatment <- factor(TAW_pupalheatexp_pupalmort$Treatment)
TAW_pupalheatexp_pupalmort$Cage_Setup <- factor(TAW_pupalheatexp_pupalmort$Cage_Setup)
TAW_pupalheatexp_pupalmort$Replicate <- factor(TAW_pupalheatexp_pupalmort$Replicate)
TAW_pupalheatexp_pupalmort$Treated <- factor(TAW_pupalheatexp_pupalmort$Treated)
TAW_pupalheatexp_pupalmort$Treated_Period <- factor(TAW_pupalheatexp_pupalmort$Treated_Period)
TAW_pupalheatexp_pupalmort$Sex <- factor(TAW_pupalheatexp_pupalmort$Sex)
str(TAW_pupalheatexp_pupalmort)## tibble [206 × 10] (S3: tbl_df/tbl/data.frame)
## $ Treatment : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 2 2 ...
## $ Period : Factor w/ 12 levels "0","1","2","3",..: 1 1 1 1 1 1 1 1 2 2 ...
## $ Cage_Setup : Factor w/ 4 levels "0","1","2","3": 1 1 1 1 1 1 1 1 2 2 ...
## $ Replicate : Factor w/ 4 levels "1","2","3","4": 1 2 3 4 1 2 3 4 1 2 ...
## $ Treated : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 2 2 ...
## $ Sex : Factor w/ 2 levels "1","2": 1 1 1 1 2 2 2 2 1 1 ...
## $ Total_Emerged : num [1:206] 5 5 5 5 8 8 8 8 5 2 ...
## $ Non_Emerging : num [1:206] 0 0 0 0 0 0 0 0 0 0 ...
## $ Percent_Nonemerging: num [1:206] 0 0 0 0 0 0 0 0 0 0 ...
## $ Treated_Period : Factor w/ 12 levels "0","1","10","11",..: 1 1 1 1 1 1 1 1 2 2 ...Subsetting so it only uses the 48h data
pup_mort_48h <- subset(TAW_pupalheatexp_pupalmort, Treatment== "0"|Treatment=="2", select=c("Period","Cage_Setup", "Sex", "Treated","Total_Emerged", "Non_Emerging", "Treated_Period", "Treatment"))If non emerging and total emerged do not sum to 8 for males and 5 for females it means there were moth escapes, these would then count as emerged. Total emerged was the number of deceased moths obtained by the end of the replicate lifespan.
Testing whether Cage_Setup influenced mortality rate of pupae.
femalefull <- subset(pup_mort_48h, Sex== "1", select=c("Period", "Cage_Setup","Treated","Total_Emerged", "Non_Emerging", "Treated_Period", "Treatment"))
femalefull$total <- 5
femalefull$rate <- femalefull$Non_Emerging/femalefull$total
malefull <- subset(pup_mort_48h, Sex== "2", select=c("Period", "Cage_Setup","Treated","Total_Emerged", "Non_Emerging", "Treated_Period","Treatment" ))
malefull$total <- 8
malefull$rate <- malefull$Non_Emerging/malefull$total
head(femalefull)## # A tibble: 6 × 9
## Period Cage_Setup Treated Total_Emerged Non_Emerging Treated_Period Treatment
## <fct> <fct> <fct> <dbl> <dbl> <fct> <fct>
## 1 0 0 0 5 0 0 0
## 2 0 0 0 5 0 0 0
## 3 0 0 0 5 0 0 0
## 4 0 0 0 5 0 0 0
## 5 6 1 1 5 0 6 2
## 6 6 1 1 4 1 6 2
## # ℹ 2 more variables: total <dbl>, rate <dbl>head(malefull)## # A tibble: 6 × 9
## Period Cage_Setup Treated Total_Emerged Non_Emerging Treated_Period Treatment
## <fct> <fct> <fct> <dbl> <dbl> <fct> <fct>
## 1 0 0 0 8 0 0 0
## 2 0 0 0 8 0 0 0
## 3 0 0 0 8 0 0 0
## 4 0 0 0 8 0 0 0
## 5 6 1 0 6 2 0 2
## 6 6 1 0 4 4 0 2
## # ℹ 2 more variables: total <dbl>, rate <dbl>glms for including Cage_Setup
Femalefull_glm <- glm(rate ~ Period*Cage_Setup, weight = total,
data = femalefull,
family = binomial(link = "logit"))
summary(Femalefull_glm)##
## Call:
## glm(formula = rate ~ Period * Cage_Setup, family = binomial(link = "logit"),
## data = femalefull, weights = total)
##
## Coefficients: (9 not defined because of singularities)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.054e+01 3.906e+03 -0.005 0.996
## Period6 1.866e+01 3.906e+03 0.005 0.996
## Period7 1.790e+01 3.906e+03 0.005 0.996
## Period8 -6.167e-10 5.967e+03 0.000 1.000
## Period9 1.866e+01 3.906e+03 0.005 0.996
## Period10 -5.780e-09 5.967e+03 0.000 1.000
## Period11 1.790e+01 3.906e+03 0.005 0.996
## Cage_Setup1 -8.059e-15 1.464e+00 0.000 1.000
## Cage_Setup2 1.253e+00 1.220e+00 1.027 0.304
## Cage_Setup3 NA NA NA NA
## Period6:Cage_Setup1 8.574e-15 1.816e+00 0.000 1.000
## Period7:Cage_Setup1 7.673e-01 1.947e+00 0.394 0.694
## Period8:Cage_Setup1 1.866e+01 4.511e+03 0.004 0.997
## Period9:Cage_Setup1 -1.866e+01 4.511e+03 -0.004 0.997
## Period10:Cage_Setup1 1.790e+01 4.511e+03 0.004 0.997
## Period11:Cage_Setup1 NA NA NA NA
## Period6:Cage_Setup2 -7.411e-02 1.538e+00 -0.048 0.962
## Period7:Cage_Setup2 -4.855e-01 1.771e+00 -0.274 0.784
## Period8:Cage_Setup2 1.664e+01 4.511e+03 0.004 0.997
## Period9:Cage_Setup2 -1.992e+01 4.511e+03 -0.004 0.996
## Period10:Cage_Setup2 -1.253e+00 6.379e+03 0.000 1.000
## Period11:Cage_Setup2 NA NA NA NA
## Period6:Cage_Setup3 NA NA NA NA
## Period7:Cage_Setup3 NA NA NA NA
## Period8:Cage_Setup3 NA NA NA NA
## Period9:Cage_Setup3 NA NA NA NA
## Period10:Cage_Setup3 NA NA NA NA
## Period11:Cage_Setup3 NA NA NA NA
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 71.627 on 57 degrees of freedom
## Residual deviance: 42.826 on 39 degrees of freedom
## AIC: 114.66
##
## Number of Fisher Scoring iterations: 18plot(Femalefull_glm)
anova(Femalefull_glm)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: rate
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 57 71.627
## Period 6 14.5059 51 57.122 0.02447 *
## Cage_Setup 2 1.7384 49 55.383 0.41928
## Period:Cage_Setup 10 12.5576 39 42.826 0.24947
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Malefull_glm <- glm(rate ~ Period*Cage_Setup, weight = total,
data = malefull,
family = binomial(link = "logit"))
summary(Malefull_glm)##
## Call:
## glm(formula = rate ~ Period * Cage_Setup, family = binomial(link = "logit"),
## data = malefull, weights = total)
##
## Coefficients: (9 not defined because of singularities)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.092e+01 3.749e+03 -0.006 0.9955
## Period6 1.959e+01 3.749e+03 0.005 0.9958
## Period7 1.898e+01 3.749e+03 0.005 0.9960
## Period8 1.026e-08 5.726e+03 0.000 1.0000
## Period9 -2.615e-08 5.726e+03 0.000 1.0000
## Period10 1.931e+01 3.749e+03 0.005 0.9959
## Period11 1.853e+01 3.749e+03 0.005 0.9961
## Cage_Setup1 1.063e+00 8.934e-01 1.190 0.2341
## Cage_Setup2 1.705e+00 8.561e-01 1.991 0.0465 *
## Cage_Setup3 NA NA NA NA
## Period6:Cage_Setup1 -6.152e-01 1.119e+00 -0.550 0.5825
## Period7:Cage_Setup1 -4.520e-01 1.197e+00 -0.378 0.7056
## Period8:Cage_Setup1 1.792e+01 4.329e+03 0.004 0.9967
## Period9:Cage_Setup1 1.746e+01 4.329e+03 0.004 0.9968
## Period10:Cage_Setup1 -1.399e+00 1.216e+00 -1.151 0.2499
## Period11:Cage_Setup1 NA NA NA NA
## Period6:Cage_Setup2 -1.468e+00 1.099e+00 -1.336 0.1815
## Period7:Cage_Setup2 -1.094e+00 1.169e+00 -0.936 0.3494
## Period8:Cage_Setup2 1.682e+01 4.329e+03 0.004 0.9969
## Period9:Cage_Setup2 1.727e+01 4.329e+03 0.004 0.9968
## Period10:Cage_Setup2 -2.102e+01 4.329e+03 -0.005 0.9961
## Period11:Cage_Setup2 NA NA NA NA
## Period6:Cage_Setup3 NA NA NA NA
## Period7:Cage_Setup3 NA NA NA NA
## Period8:Cage_Setup3 NA NA NA NA
## Period9:Cage_Setup3 NA NA NA NA
## Period10:Cage_Setup3 NA NA NA NA
## Period11:Cage_Setup3 NA NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 88.127 on 57 degrees of freedom
## Residual deviance: 39.551 on 39 degrees of freedom
## AIC: 154.46
##
## Number of Fisher Scoring iterations: 18plot(Malefull_glm)
anova(Malefull_glm)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: rate
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 57 88.127
## Period 6 27.3919 51 60.735 0.0001223 ***
## Cage_Setup 2 4.5439 49 56.192 0.1031128
## Period:Cage_Setup 10 16.6404 39 39.551 0.0827084 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Repeating the above without Cage_Setup
female_pup <- subset(pup_mort_48h, Sex== "1", select=c("Period", "Treated","Total_Emerged", "Non_Emerging", "Treated_Period", "Treatment"))
female_pup$total <- 5
female_pup$rate <- female_pup$Non_Emerging/female_pup$total
male_pup <- subset(pup_mort_48h, Sex== "2", select=c("Period", "Treated","Total_Emerged", "Non_Emerging", "Treated_Period","Treatment" ))
male_pup$total <- 8
male_pup$rate <- male_pup$Non_Emerging/male_pup$total
head(female_pup)## # A tibble: 6 × 8
## Period Treated Total_Emerged Non_Emerging Treated_Period Treatment total rate
## <fct> <fct> <dbl> <dbl> <fct> <fct> <dbl> <dbl>
## 1 0 0 5 0 0 0 5 0
## 2 0 0 5 0 0 0 5 0
## 3 0 0 5 0 0 0 5 0
## 4 0 0 5 0 0 0 5 0
## 5 6 1 5 0 6 2 5 0
## 6 6 1 4 1 6 2 5 0.2head(male_pup)## # A tibble: 6 × 8
## Period Treated Total_Emerged Non_Emerging Treated_Period Treatment total rate
## <fct> <fct> <dbl> <dbl> <fct> <fct> <dbl> <dbl>
## 1 0 0 8 0 0 0 8 0
## 2 0 0 8 0 0 0 8 0
## 3 0 0 8 0 0 0 8 0
## 4 0 0 8 0 0 0 8 0
## 5 6 0 6 2 0 2 8 0.25
## 6 6 0 4 4 0 2 8 0.5Now looking at each sex individually across all my cage types to see global pupal mortality.
FemaleMortality_glm <- glm(rate ~ Treated_Period, weight = total,
data = female_pup,
family = binomial(link = "logit"))
summary(FemaleMortality_glm)##
## Call:
## glm(formula = rate ~ Treated_Period, family = binomial(link = "logit"),
## data = female_pup, weights = total)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.197e+00 3.178e-01 -6.913 4.73e-12 ***
## Treated_Period10 -1.170e+00 1.066e+00 -1.098 0.272
## Treated_Period11 -4.418e-01 7.980e-01 -0.554 0.580
## Treated_Period6 3.254e-01 6.241e-01 0.521 0.602
## Treated_Period7 1.589e-15 6.866e-01 0.000 1.000
## Treated_Period8 -4.418e-01 7.980e-01 -0.554 0.580
## Treated_Period9 -4.418e-01 7.980e-01 -0.554 0.580
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 71.627 on 57 degrees of freedom
## Residual deviance: 68.737 on 51 degrees of freedom
## AIC: 116.57
##
## Number of Fisher Scoring iterations: 5plot(FemaleMortality_glm)
anova(FemaleMortality_glm)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: rate
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 57 71.627
## Treated_Period 6 2.8906 51 68.737 0.8225MaleMortality_glm <- glm(rate ~ Treated_Period, weight = total,
data = male_pup,
family = binomial(link = "logit"))
summary(MaleMortality_glm)##
## Call:
## glm(formula = rate ~ Treated_Period, family = binomial(link = "logit"),
## data = male_pup, weights = total)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.7984 0.2159 -8.329 <2e-16 ***
## Treated_Period10 -0.5995 0.5651 -1.061 0.2888
## Treated_Period11 0.4634 0.4159 1.114 0.2651
## Treated_Period6 0.5854 0.4057 1.443 0.1490
## Treated_Period7 0.1890 0.4434 0.426 0.6700
## Treated_Period8 -1.3371 0.7539 -1.774 0.0761 .
## Treated_Period9 -0.9096 0.6342 -1.434 0.1515
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 88.127 on 57 degrees of freedom
## Residual deviance: 74.017 on 51 degrees of freedom
## AIC: 164.93
##
## Number of Fisher Scoring iterations: 5anova(MaleMortality_glm)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: rate
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 57 88.127
## Treated_Period 6 14.11 51 74.017 0.02843 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(MaleMortality_glm)## Analysis of Deviance Table (Type II tests)
##
## Response: rate
## LR Chisq Df Pr(>Chisq)
## Treated_Period 14.11 6 0.02843 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1logit.malemort <- contrast(emmeans(MaleMortality_glm, ~ Treated_Period), adjust ="holm", method = "dunnett")
summary(logit.malemort)## contrast estimate SE df z.ratio p.value
## Treated_Period10 - Treated_Period0 -0.599 0.565 Inf -1.061 0.7954
## Treated_Period11 - Treated_Period0 0.463 0.416 Inf 1.114 0.7954
## Treated_Period6 - Treated_Period0 0.585 0.406 Inf 1.443 0.7450
## Treated_Period7 - Treated_Period0 0.189 0.443 Inf 0.426 0.7954
## Treated_Period8 - Treated_Period0 -1.337 0.754 Inf -1.774 0.4568
## Treated_Period9 - Treated_Period0 -0.910 0.634 Inf -1.434 0.7450
##
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: holm method for 6 testsCreating a summary dataset of mortality rate.
malepupsum <- summarySE(data = male_pup, measurevar = "rate", groupvars = 'Treated_Period')
malepupsum## Treated_Period N rate sd se ci
## 1 0 22 0.14204545 0.15086780 0.03216512 0.06689104
## 2 10 6 0.08333333 0.10206207 0.04166667 0.10710758
## 3 11 6 0.20833333 0.21889876 0.08936504 0.22972016
## 4 6 6 0.22916667 0.14613065 0.05965759 0.15335471
## 5 7 6 0.16666667 0.06454972 0.02635231 0.06774078
## 6 8 6 0.04166667 0.06454972 0.02635231 0.06774078
## 7 9 6 0.06250000 0.10458250 0.04269563 0.10975261femalepupsum <- summarySE(data = female_pup, measurevar = 'rate', groupvars = 'Treated_Period')
femalepupsum## Treated_Period N rate sd se ci
## 1 0 22 0.10000000 0.20236695 0.04314478 0.08972448
## 2 10 6 0.03333333 0.08164966 0.03333333 0.08568606
## 3 11 6 0.06666667 0.10327956 0.04216370 0.10838525
## 4 6 6 0.13333333 0.10327956 0.04216370 0.10838525
## 5 7 6 0.10000000 0.10954451 0.04472136 0.11495991
## 6 8 6 0.06666667 0.10327956 0.04216370 0.10838525
## 7 9 6 0.06666667 0.16329932 0.06666667 0.17137212