4.1 Additive genetics effects

In order to simulate breeding values (additive genetic effects), we can provide the simulate_population() function with the relatedness structure in the population. The simplest way to do this is providing a pedigree using the the pedigree argument (a genetic relatedness matrix could also be given to the cov_str argument). The input to this argument needs to be a list, and the name of the pedigree in the list links it with the item in the parameter list.

NOTE the simulate_population function has very little error checking of pedigree structure at the moment

When simulating breeding values, all individuals in pedigree need to be in the data_structure and vice versa. Having unsampled individuals (for example the base population) can be achieved in the sampling stage (not implemented yet).

Lets start by importing a pedigree

library(MCMCglmm)
data(BTped)
head(BTped)

##    animal  dam sire
## 1 R187557 <NA> <NA>
## 2 R187559 <NA> <NA>
## 3 R187568 <NA> <NA>
## 4 R187518 <NA> <NA>
## 5 R187528 <NA> <NA>
## 6 R187945 <NA> <NA>

We can use this pedigree as a data_structure

squid_data <- simulate_population(
  data_structure = BTped, 
 
  pedigree = list(animal=BTped),
  
  parameters =list(
    animal = list(
      vcov = 0.2
    ),
    residual = list(
      vcov = 0.5
    )
  )
)

data <- get_population_data(squid_data)
head(data)

##                  y animal_effect   residual  animal  dam sire squid_pop
## R187557 -0.7104031    -0.1749195 -0.5354836 R187557 <NA> <NA>         1
## R187559  1.2892857     0.6665992  0.6226864 R187559 <NA> <NA>         1
## R187568  1.6393384     0.6900526  0.9492857 R187568 <NA> <NA>         1
## R187518  1.2545103     0.3999815  0.8545287 R187518 <NA> <NA>         1
## R187528 -0.4945722    -0.3827869 -0.1117853 R187528 <NA> <NA>         1
## R187945  0.1780773    -0.3703026  0.5483799 R187945 <NA> <NA>         1

# Ainv<-inverseA(BTped)$Ainv
# mod <- MCMCglmm(y~1, random=~ animal,data=data,ginverse=list(animal=Ainv),verbose=FALSE)
# summary(mod)

We might want to simulate repeated measurements to allow estimation of permanent environment effects. This is where being able to have something in the parameter list with a different name to the grouping factor is useful. In this way permanent environmental and additive genetic effects can be simulated in different parts of the parameter list, and linked to the same part of the data_structure.

## make data structure with two observations per individual
ds <- data.frame(individual=rep(BTped[,1], 2))

squid_data <- simulate_population(
  data_structure = ds, 
  pedigree=list(animal=BTped),
  parameters = list(
    individual = list(
      vcov = 0.3
    ),
    animal = list(
      group="individual",
      vcov = 0.2
    ),
    residual = list(
      vcov = 0.5
    )
  )
)

data <- get_population_data(squid_data)
head(data)

##                   y individual_effect animal_effect   residual individual
## R187888 -0.29800108        -0.3030617    -0.2605827  0.2656433    R187557
## R187646  0.17464214         0.3344371    -0.5309132  0.3711182    R187559
## R187330  1.17443051         0.7872808     0.9885985 -0.6014488    R187568
## R187374 -0.88713844        -0.4727921    -0.1421482 -0.2721982    R187518
## R187225 -0.06469935        -0.2399442     0.5130975 -0.3378527    R187528
## R187133  1.72913172         0.5764696     0.7417825  0.4108797    R187945
##         squid_pop
## R187888         1
## R187646         1
## R187330         1
## R187374         1
## R187225         1
## R187133         1

# Ainv<-inverseA(BTped)$Ainv
# data$animal_id <- data$individual
# mod <- MCMCglmm(y~1, random=~ individual + animal_id,data=data,ginverse=list(animal_id=Ainv),verbose=FALSE)
# summary(mod)