In this vignette, the `fullfact`

package is explored using
advanced (functions designated by the number 2) for the standard model,
i.e. containing random effects for dam, sire, and dam by sire, for
normal data type or error structure, with the options of including
additional random effects for one position (e.g. tank) and/or one block
effect (e.g. several blocks of 2 \(\times\) 2 factorial matings).

Simple (functions designated by no number) for the standard model
only is explored in the vignette **Simple Normal Data
Example**.

Expert (functions designated by the number 3) for the standard model
with the ability of the user to include additional fixed and/or random
effects, such as a model including environment treatments and their
interactions is explored in the vignette **Expert Normal Data
Example**.

Non-normal error structures (e.g. binary, proportion, and/or count data types) are explored in another three vignettes: (1) Simple Non-normal Data Example, (2) Advanced Non-normal Data Example, and (3) Expert Non-normal Data Example.

The example data set is an 11 \(\times\) 11 full factorial mating: 11 dams and 11 sires with all combinations resulting in 121 families. There are 10 observations per family or 5 observations for each of two replicates per family.

```
library("fullfact")
data(chinook_length)
head(chinook_length)
#> family repli dam sire tray cell length egg_size
#> 1 f1 r1 d1 s1 t7 1A 22.7 7.27
#> 2 f1 r2 d1 s1 t8 1A 22.0 7.27
#> 3 f1 r1 d1 s1 t7 1A 23.8 7.27
#> 4 f1 r2 d1 s1 t8 1A 21.8 7.27
#> 5 f1 r1 d1 s1 t7 1A 22.6 7.27
#> 6 f1 r2 d1 s1 t8 1A 22.4 7.27
```

Displayed are columns for family identities (ID), replicate ID, dam ID, sire ID, incubation tray ID, incubation cell ID (within tray), Chinook salmon length (mm) at hatch, and dam egg size (mm).

Model random effects are dam, sire, and dam by sire. Options to include one random position and/or one random block effect(s). Extracts the dam, sire, dam, dam by sire, and residual variance components. Calculates the total variance component. Calculates the additive genetic, non-additive genetic, and maternal variance components. Extracts optional position and block variance components.

Assuming the effects of epistasis are of negligible importance, the
additive genetic variance (V_{A}) component is calculated as
four times the sire (V_{S}), the non-additive genetic variance
(V_{N}) component as four times the dam by sire interaction
(V_{D\(\times\)S}), and the
maternal variance component (V_{M}) as the dam (V_{D}) –
sire (V_{S}) (Lynch and Walsh 1998, p. 603). When there is
epistasis, those variance components will be overestimated and this may
explain why the percentage of phenotypic variance explained by the
components can add up to more than 100% in certain cases.

`position`

is the column name containing position factor
information.

`block`

is the column name containing block factor
information.

Significance values for the random effects are determined using likelihood ratio tests (Bolker et al. 2009).

Lynch M, Walsh B. 1998. Genetics and Analysis of Quantitative Traits. Sinauer Associates, Massachusetts.

Bolker BM, Brooks ME, Clark CJ, Geange SW, Poulsen JR, Stevens MHH, White J-SS. 2009. Generalized linear mixed models: a practical guide for ecology and evolution. Trends in Ecology and Evolution 24(3): 127-135. DOI: 10.1016/j.tree.2008.10.008

For this example, we explore position (tray) effects for Chinook salmon length.

```
length_mod2<- observLmer2(observ=chinook_length,dam="dam",sire="sire",response="length",
position="tray")
#> [1] "2024-02-04 15:48:48 PST"
#> Time difference of 0.3131449 secs
length_mod2
#> $random
#> effect variance percent d.AIC d.BIC Chi.sq p.value
#> 1 dam:sire 0.1798015 17.28121 141.49865 136.400272 143.49865 4.573043e-33
#> 2 tray 0.1077362 10.35482 121.50329 116.404911 123.50329 1.082069e-28
#> 3 sire 0.0000000 0.00000 -2.00000 -7.098376 0.00000 1.000000e+00
#> 4 dam 0.2029810 19.50905 29.98924 24.890867 31.98924 1.550287e-08
#>
#> $other
#> component variance percent
#> 1 Residual 0.5499266 52.85493
#> 2 Total 1.0404453 100.00000
#>
#> $calculation
#> component variance percent
#> 1 additive 0.0000000 0.00000
#> 2 nonadd 0.7192061 69.12483
#> 3 maternal 0.2029810 19.50905
```

Produces a list object containing three data frames. Each data frame contains the raw variance components and the variance components as a percentage of the total variance component. The first data frame also contains the difference in AIC and BIC, and likelihood ratio test Chi-square and p-value for all random effects.

*Note*

Default is Restricted maximum likelihood (REML) as
`ml = F`

. Option for maximum likelihood (ML) is
`ml = T`

.

Maximum likelihood (ML) estimates the parameters that maximize the likelihood of the observed data and has the advantage of using all the data and accounting for non-independence (Lynch and Walsh 1998, p. 779; Bolker et al. 2009). On the other hand, ML has the disadvantage of assuming that all fixed effects are known without error, producing a downward bias in the estimation of the residual variance component. This bias can be large if there are lots of fixed effects, especially if sample sizes are small. Restricted maximum likelihood (REML) has the advantage of not assuming the fixed effects are known and averages over the uncertainty, so there can be less bias in the estimation of the residual variance component. However, REML only maximizes a portion of the likelihood to estimate the effect parameters, but is the preferred method for analyzing large data sets with complex structure.

Power values are calculated by stochastically simulating data for a number of iterations and then calculating the proportion of P-values less than \(\alpha\) (e.g. 0.05) for each component (Bolker 2008). Simulated data are specified by inputs for variance component values and the sample sizes.

Bolker BM. 2008. Ecological Models and Data in R. Princeton University Press, Princeton.

Defaults are `alpha = 0.05`

for 5%,
`nsim = 100`

for 100 simulations, and `ml = F`

for
REML.

`varcomp`

is a vector of dam, sire, dam by sire, residual,
position and/or block variance components, i.e. c(dam,sire,dam \(\times\) sire,residual,position/block). If
there is a position and a block, c(dam,sire,dam \(\times\) sire,residual,position,block).

`nval`

is a vector of dam, sire, offspring per family, and
offspring per position or number of block sample sizes,
i.e. c(dam,sire,offspring,position/block). If there is a position and a
block, c(dam,sire,offspring,position,block).

`position`

is optional number of positions.

`block`

is optional vector of dams and sires per block,
e.g. c(2,2).

For this example, the variance components of `observLmer2`

above are used (i.e. dam= 0.2030, sire= 0, dam \(\times\) sire= 0.1798, residual= 0.5499,
tray= 0.1077) and the sample size of the Chinook salmon data set
(i.e. dam= 11, sire= 11, offspring= 10, offspring per position= 110).
Position was represented by 11 trays. The actual design was composed of
16 trays with 55–80 offspring each. However, `powerLmer2`

uses an equal number of offspring per position, so the number of trays
was decreased from 16 to 11.

Full analysis is 100 simulations. Example has 25 simulations.

```
#>powerLmer2(varcomp=c(0.2030,0,0.1798,0.5499,0.1077),nval=c(11,11,10,110),position=11) #full
#25 simulations
powerLmer2(varcomp=c(0.2030,0,0.1798,0.5499,0.1077),nval=c(11,11,10,110),position=11,nsim=25)
#> [1] "2024-02-04 15:48:49 PST"
#> [1] "Starting simulation: 1"
#> [1] "Starting simulation: 2"
#> [1] "Starting simulation: 3"
#> [1] "Starting simulation: 4"
#> [1] "Starting simulation: 5"
#> [1] "Starting simulation: 6"
#> [1] "Starting simulation: 7"
#> [1] "Starting simulation: 8"
#> [1] "Starting simulation: 9"
#> [1] "Starting simulation: 10"
#> [1] "Starting simulation: 11"
#> [1] "Starting simulation: 12"
#> [1] "Starting simulation: 13"
#> [1] "Starting simulation: 14"
#> [1] "Starting simulation: 15"
#> [1] "Starting simulation: 16"
#> [1] "Starting simulation: 17"
#> [1] "Starting simulation: 18"
#> [1] "Starting simulation: 19"
#> [1] "Starting simulation: 20"
#> [1] "Starting simulation: 21"
#> [1] "Starting simulation: 22"
#> [1] "Starting simulation: 23"
#> [1] "Starting simulation: 24"
#> [1] "Starting simulation: 25"
#> Time difference of 8.156208 secs
#> term n var_in var_out power
#> 1 dam 11 0.2030 0.17589809 1.00
#> 2 sire 11 0.0000 0.00631872 0.08
#> 3 dam.sire 121 0.1798 0.16449326 1.00
#> 4 position 11 0.1077 0.09609158 1.00
#> 5 residual NA 0.5499 0.54673647 NA
#Block examples using 8 dams, 8 sires (as four 2x2 blocks), and 20 offspring per family
#>powerLmer2(varcomp=c(0.2030,0,0.1798,0.5499,0.1077),nval=c(8,8,20,4),block=c(2,2))
#>powerLmer2(varcomp=c(0.2030,0,0.1798,0.5499,0.1077,0.1077),nval=c(8,8,20,40,4),position=8,
#>block=c(2,2)) #with position
```

There is sufficient power (\(\ge\) 0.8) for dam and dam by sire variance components, whereas there is insufficient power (< 0.8) for the sire variance component. Albeit, the sire component is near zero, so the low power may be an artifact. There was also sufficient power for position (tray) variance components. In the cases of insufficient power, the sample size of dam, sire, and/or offspring can be increased until there is sufficient power.

Taking the reverse approach (can the sample size of dam, sire, or offspring be reduced while maintaining sufficient power?) using the same variance components and offspring sample size, dam and sire sample sizes could be reduced from 11 to 7. The position number and its sample size was reduced accordingly, i.e. 7 dams \(\times\) 7 sires \(\times\) 10 offspring = 490, divided by 10 trays for 49 offspring each.

```
#>powerLmer2(varcomp=c(0.2030,0,0.1798,0.5499,0.1077),nval=c(7,7,10,49),position=10) #full
#25 simulations
powerLmer2(varcomp=c(0.2030,0,0.1798,0.5499,0.1077),nval=c(7,7,10,49),position=10,nsim=25)
#> [1] "2024-02-04 15:48:57 PST"
#> [1] "Starting simulation: 1"
#> [1] "Starting simulation: 2"
#> [1] "Starting simulation: 3"
#> [1] "Starting simulation: 4"
#> [1] "Starting simulation: 5"
#> [1] "Starting simulation: 6"
#> [1] "Starting simulation: 7"
#> [1] "Starting simulation: 8"
#> [1] "Starting simulation: 9"
#> [1] "Starting simulation: 10"
#> [1] "Starting simulation: 11"
#> [1] "Starting simulation: 12"
#> [1] "Starting simulation: 13"
#> [1] "Starting simulation: 14"
#> [1] "Starting simulation: 15"
#> [1] "Starting simulation: 16"
#> [1] "Starting simulation: 17"
#> [1] "Starting simulation: 18"
#> [1] "Starting simulation: 19"
#> [1] "Starting simulation: 20"
#> [1] "Starting simulation: 21"
#> [1] "Starting simulation: 22"
#> [1] "Starting simulation: 23"
#> [1] "Starting simulation: 24"
#> [1] "Starting simulation: 25"
#> Time difference of 5.511813 secs
#> term n var_in var_out power
#> 1 dam 7 0.2030 0.22885061 0.84
#> 2 sire 7 0.0000 0.00587355 0.00
#> 3 dam.sire 49 0.1798 0.17430060 1.00
#> 4 position 10 0.1077 0.09689070 1.00
#> 5 residual NA 0.5499 0.55339969 NA
```

Confidence intervals for the additive genetic, non-additive genetic, and maternal variance components can be produced using the bootstrap-t resampling method described by Efron and Tibshirani (1993, p. 160‒162). Observations are resampled with replacement until the original sample size is reproduced. The resampled data are then used in the model and the additive genetic, non-additive genetic, and maternal variance components are extracted. The process is repeated for a number of iterations, typically 1,000 times, to produce a distribution for each component. The confidence interval lower and upper limits and median are extracted from the distribution.

Efron B, Tibshirani R. 1993. An Introduction to the Bootstrap. Chapman and Hall, New York.

The `resampRepli`

function is used to bootstrap resample
observations grouped by replicate ID within family ID for a specified
number of iterations to create the resampled data set. A similar
`resampFamily`

function is able to resample observations
grouped by family ID only.

`copy`

is a vector of column numbers (to copy the
contents). Does not need to contain the family and/or replicate
columns.

Full analysis is 1000 iterations. Example has 5 iterations.

```
#>resampRepli(dat=chinook_length,copy=c(3:8),family="family",replicate="repli",iter=1000) #full
#>resampFamily(dat=chinook_length,copy=c(3:8),family="family",iter=1000) #family only
resampRepli(dat=chinook_length,copy=c(3:8),family="family",replicate="repli",iter=5) #5 iterations
```

Because of the large file sizes that can be produced, the resampling of each replicate Y per family X is saved separately as a common separated (X_Y_resampR.csv) file in the working directory. These files are merged to create the final resampled data set (resamp_datR.csv).

If using `resampFamily`

, the file names are X_resampF.csv
per family and resamp_datF.csv for the final resampled data set.

The equivalent to `observLmer2`

is available for the final
bootstrap resampled data set, i.e. `resampLmer2`

.

Default is `ml = F`

for REML. The starting model number
`start =`

and ending model number `end =`

need to
be specified.

Full analysis is 1000 iterations. Example has 5 iterations.

```
#>length_datR<- read.csv("resamp_datR.csv") #1000 iterations
#>length_rcomp2<- resampLmer2(resamp=length_datR,dam="dam",sire="sire",response="length",
#>position="tray",start=1,end=1000) #full
data(chinook_resampL) #5 iterations
head(chinook_resampL)
#> dam1 sire1 tray1 cell1 length1 egg_size1 dam2 sire2 tray2 cell2 length2
#> 1 d1 s1 t7 1A 23.1 7.27 d1 s1 t7 1A 23.1
#> 2 d1 s1 t7 1A 22.7 7.27 d1 s1 t7 1A 22.7
#> 3 d1 s1 t7 1A 22.9 7.27 d1 s1 t7 1A 22.9
#> 4 d1 s1 t7 1A 23.1 7.27 d1 s1 t7 1A 22.7
#> 5 d1 s1 t7 1A 23.1 7.27 d1 s1 t7 1A 22.9
#> 6 d1 s1 t8 1A 22.4 7.27 d1 s1 t8 1A 22.3
#> egg_size2 dam3 sire3 tray3 cell3 length3 egg_size3 dam4 sire4 tray4 cell4
#> 1 7.27 d1 s1 t7 1A 23.8 7.27 d1 s1 t7 1A
#> 2 7.27 d1 s1 t7 1A 23.8 7.27 d1 s1 t7 1A
#> 3 7.27 d1 s1 t7 1A 23.8 7.27 d1 s1 t7 1A
#> 4 7.27 d1 s1 t7 1A 23.8 7.27 d1 s1 t7 1A
#> 5 7.27 d1 s1 t7 1A 23.8 7.27 d1 s1 t7 1A
#> 6 7.27 d1 s1 t8 1A 22.4 7.27 d1 s1 t8 1A
#> length4 egg_size4 dam5 sire5 tray5 cell5 length5 egg_size5
#> 1 22.6 7.27 d1 s1 t7 1A 22.9 7.27
#> 2 23.8 7.27 d1 s1 t7 1A 22.6 7.27
#> 3 23.1 7.27 d1 s1 t7 1A 22.6 7.27
#> 4 22.9 7.27 d1 s1 t7 1A 22.9 7.27
#> 5 23.8 7.27 d1 s1 t7 1A 22.7 7.27
#> 6 22.0 7.27 d1 s1 t8 1A 22.4 7.27
length_rcomp2<- resampLmer2(resamp=chinook_resampL,dam="dam",sire="sire",response="length",
position="tray",start=1,end=5)
#> [1] "2024-02-04 15:49:03 PST"
#> [1] "Working on model: 1"
#> [1] "Working on model: 2"
#> [1] "Working on model: 3"
#> [1] "Working on model: 4"
#> [1] "Working on model: 5"
#> Time difference of 0.4106691 secs
length_rcomp2[1:5,]
#> dam:sire tray sire dam Residual Total additive
#> 1 0.2156343 0.1010264 9.658700e-09 0.1933548 0.4996255 1.0096410 3.863480e-08
#> 2 0.2423708 0.1084802 3.551509e-10 0.1954621 0.5204227 1.0667358 1.420604e-09
#> 3 0.2204124 0.1072239 3.857906e-10 0.1713985 0.4844660 0.9835008 1.543162e-09
#> 4 0.1959125 0.1204563 1.996910e-03 0.1931546 0.4919967 1.0035171 7.987638e-03
#> 5 0.2263470 0.1308195 0.000000e+00 0.2085843 0.4787757 1.0445264 0.000000e+00
#> nonadd maternal
#> 1 0.8625374 0.1933548
#> 2 0.9694834 0.1954621
#> 3 0.8816494 0.1713985
#> 4 0.7836500 0.1911577
#> 5 0.9053880 0.2085843
```

The function provides a data frame with columns containing the raw variance components for dam, sire, dam by sire, residual, total, additive genetic, non-additive genetic, and maternal. Also columns containing the raw variance components for the options of position and/or block. The number of rows in the data frame matches the number of iterations in the resampled data set and each row represents a model number.

Extract the bootstrap-t confidence intervals (CI) and median for the
additive genetic, non-additive genetic, and maternal values from the
data frame of models produced using `resampLmer2`

. Also
extracts intervals for optional position and block variance
components.

Default confidence interval is 95% as `level = 95`

.

```
#>ciMANA(comp=length_rcomp1) #full
data(chinook_bootL) #same as length_rcomp2 1000 models
ciMANA2(comp=chinook_bootL,position="tray")
#> $raw
#> component lower median upper
#> 1 additive 0.000 0.000 0.030
#> 2 nonadd 0.713 0.844 0.996
#> 3 maternal 0.168 0.201 0.233
#> 4 tray 0.091 0.112 0.133
#>
#> $percentage
#> component lower median upper
#> 1 additive 0.0 0.0 2.9
#> 2 nonadd 69.4 81.3 94.2
#> 3 maternal 16.7 19.3 22.0
#> 4 tray 9.0 10.8 12.6
```

The raw values are presented and are converted to a percentage of the
total variance for each model. Defaults are the number of decimal places
to round CI raw values as `rnd_r = 3`

and to round the CI
percent values as `rnd_p = 1`

.

The bootstrap-t method may produce medians that are largely different from the observed values. In the present example, the raw 95% CI are a good fit to the observed values. Nonetheless, options are provided below for 95% CI that are a poor fit.

The BCa method (bias and acceleration) described by Efron and Tibshirani (1993, p.184‒188) can be used for the correction of bootstrap-t confidence intervals.

`bias`

is a vector of additive, non-additive, maternal,
position and/or block variance components,
i.e. c(additive,non-additive,maternal,position/block), from the raw
observed variance components of `observLmer2`

. If there is a
position and a block,
c(additive,non-additive,maternal,position,block).

For this example, additive= 0, non-additive= 0.7192, maternal=
0.2030, and position= 0.1077 from `observLmer2`

above.

The ‘bias fail’ warning is if the bias calculation is infinity (negative or positive), e.g. bias contains a zero value, so the uncorrected confidence interval is displayed for the component.

```
ciMANA2(comp=chinook_bootL,position="tray",bias=c(0,0.7192,0.2030,0.1077)) #bias only
#> $raw
#> component lower median upper change
#> 1 additive 0.000 0.000 0.030 bias fail
#> 2 nonadd 0.619 0.624 0.735 <NA>
#> 3 maternal 0.174 0.206 0.240 <NA>
#> 4 tray 0.081 0.103 0.125 <NA>
#>
#> $percentage
#> component lower median upper change
#> 1 additive 0.0 0.0 2.9 bias fail
#> 2 nonadd 62.5 62.6 71.5 <NA>
#> 3 maternal 17.0 19.7 22.5 <NA>
#> 4 tray 8.0 10.0 11.9 <NA>
```

`accel`

for acceleration correction uses the delete-one
observation jackknife data set. See next section.

```
data(chinook_jackL)
#bias and acceleration
ciMANA2(comp=chinook_bootL,position="tray",bias=c(0,0.7192,0.2030,0.1077),accel=chinook_jackL)
#> $raw
#> component lower median upper change
#> 1 additive 0.000 0.000 0.030 bias fail
#> 2 nonadd 0.619 0.624 0.735 <NA>
#> 3 maternal 0.174 0.206 0.240 <NA>
#> 4 tray 0.081 0.103 0.125 <NA>
#>
#> $percentage
#> component lower median upper change
#> 1 additive 0.0 0.0 2.9 bias fail
#> 2 nonadd 62.5 62.6 71.5 <NA>
#> 3 maternal 17.0 19.7 22.5 <NA>
#> 4 tray 8.0 10.0 11.9 <NA>
```

Jackknife resampling is another method for producing confidence intervals.

The equivalent to `observLmer2`

is available for jackknife
resampling, i.e. `JackLmer2`

, using the observed data
frame.

Default is delete-one jackknife resampling as `size = 1`

and REML as `ml = F`

.

Full analysis uses all observations. Example has the first 10 observations.

```
#>length_jack2<- JackLmer2(observ=chinook_length,dam="dam",sire="sire",response="length",
#>position="tray") #full, all observations
length_jack2<- JackLmer2(observ=chinook_length,dam="dam",sire="sire",response="length",
position="tray",first=10) #first 10 observations
#> [1] "2024-02-04 15:49:03 PST"
#> [1] "Removing observation: 1 of 1210"
#> [1] "Removing observation: 2 of 1210"
#> [1] "Removing observation: 3 of 1210"
#> [1] "Removing observation: 4 of 1210"
#> [1] "Removing observation: 5 of 1210"
#> [1] "Removing observation: 6 of 1210"
#> [1] "Removing observation: 7 of 1210"
#> [1] "Removing observation: 8 of 1210"
#> [1] "Removing observation: 9 of 1210"
#> [1] "Removing observation: 10 of 1210"
#> Time difference of 0.8632441 secs
head(length_jack2)
#> dam:sire tray sire dam Residual Total additive
#> 1 0.1798331 0.1078537 0.000000e+00 0.2029751 0.5503930 1.041055 0.000000e+00
#> 2 0.1798645 0.1075831 0.000000e+00 0.2029644 0.5503499 1.040762 0.000000e+00
#> 3 0.1798327 0.1069848 6.634044e-11 0.2028003 0.5495369 1.039155 2.653617e-10
#> 4 0.1799274 0.1074530 0.000000e+00 0.2029821 0.5501989 1.040561 0.000000e+00
#> 5 0.1798594 0.1079313 8.104640e-11 0.2030077 0.5503460 1.041144 3.241856e-10
#> 6 0.1797887 0.1078271 0.000000e+00 0.2029536 0.5504038 1.040973 0.000000e+00
#> nonadd maternal
#> 1 0.7193325 0.2029751
#> 2 0.7194581 0.2029644
#> 3 0.7193308 0.2028003
#> 4 0.7197097 0.2029821
#> 5 0.7194376 0.2030077
#> 6 0.7191550 0.2029536
```

Because the delete-one observation jackknife resampling may be
computationally intensive for large data sets, the
`JackLmer2`

function has the option of delete-d observation
jackknife resampling, for which d > 1. The rows of the observed data
frame are shuffled and a block of observations of size d is deleted
sequentially. For example, delete-5 observation jackknife resampling is
specified as `size = 5`

, which deletes a block of 5
observations.

Full analysis uses all observations. Example has the first 10 observations.

```
#>length_jack2D<- JackLmer2(observ=chinook_length,dam="dam",sire="sire",response="length",
#>position="tray",size=5) #full
length_jack2D<- JackLmer2(observ=chinook_length,dam="dam",sire="sire",response="length",
position="tray",size=5,first=10) #first 10
#> [1] "2024-02-04 15:49:04 PST"
#> [1] "Removing block: 1 of 242"
#> [1] "Removing block: 2 of 242"
#> [1] "Removing block: 3 of 242"
#> [1] "Removing block: 4 of 242"
#> [1] "Removing block: 5 of 242"
#> [1] "Removing block: 6 of 242"
#> [1] "Removing block: 7 of 242"
#> [1] "Removing block: 8 of 242"
#> [1] "Removing block: 9 of 242"
#> [1] "Removing block: 10 of 242"
#> Time difference of 0.9214559 secs
head(length_jack2D)
#> dam:sire tray sire dam Residual Total additive
#> 1 0.1785926 0.1094573 3.475110e-09 0.2038849 0.5484760 1.040411 1.390044e-08
#> 2 0.1804442 0.1071383 2.540238e-10 0.2027912 0.5518833 1.042257 1.016095e-09
#> 3 0.1787291 0.1081710 0.000000e+00 0.2025798 0.5518801 1.041360 0.000000e+00
#> 4 0.1822716 0.1079838 0.000000e+00 0.2029429 0.5503571 1.043555 0.000000e+00
#> 5 0.1829877 0.1081524 0.000000e+00 0.2038243 0.5502178 1.045182 0.000000e+00
#> 6 0.1792990 0.1074311 1.546598e-10 0.2026566 0.5519269 1.041314 6.186392e-10
#> nonadd maternal
#> 1 0.7143705 0.2038849
#> 2 0.7217768 0.2027912
#> 3 0.7149162 0.2025798
#> 4 0.7290865 0.2029429
#> 5 0.7319509 0.2038243
#> 6 0.7171962 0.2026566
```

Extract the jackknife confidence intervals (CI) and median for the
additive genetic, non-additive genetic, and maternal values from the
data frame of models produced using `JackLmer2`

. Also
extracts intervals for optional position and block variance
components.

The mean and the standard error of pseudo-values for each variance component are calculated (Efron and Tibshirani 1993, p.184‒188). The standard error is then used with the Student’s t distribution to provide the lower and upper limits for the confidence interval. For delete-d jackknife resampling, M degrees of freedom were used for producing the confidence interval (Martin et al. 2004): M = N / d, where N is the total number of observations and d is the number of deleted observations. Large values of M, such as 1,000, can translate to the delete-d jackknife resampling method approaching bootstrap resampling expectations (Efron and Tibshirani 1993, p. 149).

Martin, H., Westad, F. & Martens, H. (2004). Improved Jackknife Variance Estimates of Bilinear Model Parameters. COMPSTAT 2004 – Proceedings in Computational Statistics 16th Symposium Held in Prague, Czech Republic, 2004 (ed J. Antoch), pp. 261-275. Physica-Verlag HD, Heidelberg.

Default confidence interval is 95% as `level = 95`

.

`full`

is a vector of additive, non-additive, maternal,
total, position and/or block variance components,
i.e. c(additive,non-additive,maternal,total,position/block), from the
raw observed variance components of `observLmer2`

. If there
is a position and a block,
c(additive,non-additive,maternal,total,position,block).

For this example, additive= 0, non-additive= 0.7192, maternal=
0.2030, total= 1.0404, and position= 0.1077 from
`observLmer2`

above.

```
#full, all observations
#>ciJack2(comp=length_jack2,position="tray",full=c(0,0.7192,0.2030,1.0404,0.1077))
data(chinook_jackL) #same as length_jack2 all observations
ciJack2(comp=chinook_jackL,position="tray",full=c(0,0.7192,0.2030,1.0404,0.1077))
#> $raw
#> component lower mean upper
#> 1 additive 0.000 0.000 0.000
#> 2 nonadd 0.481 0.691 0.900
#> 3 maternal 0.197 0.244 0.291
#> 4 tray 0.024 0.056 0.088
#>
#> $percentage
#> component lower mean upper
#> 1 additive 0.0 0.0 0.0
#> 2 nonadd 52.2 69.7 87.2
#> 3 maternal 20.6 24.4 28.2
#> 4 tray 3.2 5.9 8.6
```

The raw values are presented and are converted to a percentage of the
total variance for each model. Defaults are the number of decimal places
to round CI raw values as `rnd_r = 3`

and to round the CI
percent values as `rnd_p = 1`

.

The `barMANA`

and `boxMANA`

functions are
simple plotting functions for the confidence intervals or all values,
respectively, from the bootstrap and jackknife approaches. Default is to
display the percentage values as `type = perc`

. Raw values
can be displayed as `type = raw`

.

Within the functions, there are simple plot modifications available.
For the y-axis, min and max values can be species as `ymin`

and `ymax`

, as well as the increment as `yunit`

.
Also, magnification of the axis unit as `cex_yaxis`

and label
as `cex_ylab`

. The position of the legend can be specified as
`leg`

. Default is “topright”.

The `barMANA`

function produces bar graphs with the
bootstrap-t median (`ciMANA2`

) or jackknife pseudo-value mean
(`ciJack2`

) as the top of the shaded bar and error bars
covering the range of the confidence interval for each of the additive
genetic, non-additive genetic, and maternal values of a phenotypic
trait.

The length of the error bar can be specified in inches as
`bar_len`

.

```
length_ci<- ciJack2(comp=chinook_jackL,position="tray",full=c(0,0.7192,0.2030,1.0404,0.1077))
oldpar<- par(mfrow=c(2,1))
barMANA(ci_dat=length_ci) #basic, top
barMANA(ci_dat=length_ci,bar_len=0.3,yunit=20,ymax=100,cex_ylab=1.3) #modified, bottom
```

Different traits can also be combined on the same bar plot using
`trait`

specified in `ciMANA2`

or
`ciJack2`

. The information is combined into a list object.
For the example, the jackknife CI is duplicated to simulate ‘different
traits’.

```
length_ci1<- ciJack2(comp=chinook_jackL,position="tray",full=c(0,0.7192,0.2030,1.0404,0.1077),
trait="length_1")
length_ci2<- ciJack2(comp=chinook_jackL,position="tray",full=c(0,0.7192,0.2030,1.0404,0.1077),
trait="length_2")
comb_bar<- list(raw=rbind(length_ci1$raw,length_ci2$raw),
percentage=rbind(length_ci1$percentage,length_ci2$percentage))
barMANA(ci_dat=comb_bar,bar_len=0.3,yunit=20,ymax=100,cex_ylab=1.3)
```

The legend is slightly off in the presented html version but is fine with the R plotting device.

The `boxMANA`

function produces box plots using all values
for the bootstrap-t resampling data set (`resampLmer2`

) or
jackknife resampling data set (`JackLmer2`

).

```
oldpar<- par(mfrow=c(2,1))
boxMANA(comp=chinook_bootL) #from resampLmer2, basic, top
boxMANA(comp=chinook_bootL,yunit=20,ymax=100,cex_ylab=1.3,leg="topleft") #modified, bottom
```

Different traits can also be combined on the same box plot by adding a “trait” column to the resampling data set. For the example, the bootstrap-t data frame is duplicated to simulate ‘different traits’.

```
chinook_bootL1<- chinook_bootL; chinook_bootL2<- chinook_bootL #from resampLmer2
chinook_bootL1$trait<- "length_1"; chinook_bootL2$trait<- "length_2"
comb_boot<- rbind(chinook_bootL1,chinook_bootL2)
comb_boot$trait<- as.factor(comb_boot$trait)
boxMANA(comb_boot,yunit=20,ymax=100,cex_ylab=1.3,leg="topleft")
```

The recommended follow-up vignette is the **Expert Normal Data
Example**, covering the standard model with the ability of the
user to include additional fixed and/or random effects, such as a model
including environment treatments and their interactions.