The AA-Aa-aa Blood Group System

Estimates the parameter of the AA-Aa-aa blood group system, with or without Hardy Weinberg equilibrium.

AA.Aa.aa(linkp = "logitlink", linkf = "logitlink", inbreeding = FALSE,
         ipA = NULL, ifp = NULL, zero = NULL)

Arguments

linkp, linkf

Link functions applied to pA and f. See Links for more choices.

ipA, ifp

Optional initial values for pA and f.

inbreeding

Logical. Is there inbreeding?

zero

See CommonVGAMffArguments for information.

Details

This one or two parameter model involves a probability called pA. The probability of getting a count in the first column of the input (an AA) is pA*pA. When inbreeding = TRUE, an additional parameter f is used. If inbreeding = FALSE then \(f = 0\) and Hardy-Weinberg Equilibrium (HWE) is assumed. The EIM is used if inbreeding = FALSE.

Value

An object of class "vglmff" (see vglmff-class). The object is used by modelling functions such as vglm and vgam.

References

Weir, B. S. (1996). Genetic Data Analysis II: Methods for Discrete Population Genetic Data, Sunderland, MA: Sinauer Associates, Inc.

Author

T. W. Yee

Note

The input can be a 3-column matrix of counts, where the columns are AA, Ab and aa (in order). Alternatively, the input can be a 3-column matrix of proportions (so each row adds to 1) and the weights argument is used to specify the total number of counts for each row.

Warning

Setting inbreeding = FALSE makes estimation difficult with non-intercept-only models. Currently, this code seems to work with intercept-only models.

Examples

y <- cbind(53, 95, 38)
fit1 <- vglm(y ~ 1, AA.Aa.aa, trace = TRUE)
#> Iteration 1: deviance = 0.1478917
#> Iteration 2: deviance = 0.1478917
#> Iteration 3: deviance = 0.1478917
#> Taking a modified step...
#> Iteration  3 :  deviance = 0.1478917
#> Warning: some quantities such as z, residuals, SEs may be inaccurate due to convergence at a half-step
fit2 <- vglm(y ~ 1, AA.Aa.aa(inbreeding = TRUE), trace = TRUE)
#> Iteration 1: deviance = 3.430586
#> Iteration 2: deviance = 0.5959085
#> Iteration 3: deviance = 0.189806
#> Iteration 4: deviance = 0.1478993
#> Iteration 5: deviance = 0.1478917
#> Warning: 1 diagonal elements of the working weights variable 'wz' have been replaced by 1.819e-12
#> Iteration 6: deviance = 0.1478917
#> Taking a modified step....................
#> Warning: iterations terminated because half-step sizes are very small
#> Warning: some quantities such as z, residuals, SEs may be inaccurate due to convergence at a half-step
rbind(y, sum(y) * fitted(fit1))
#>         AA       Aa       aa
#>   53.00000 95.00000 38.00000
#> 1 54.30242 92.39516 39.30242
Coef(fit1)  # Estimated pA
#>        pA 
#> 0.5403226 
Coef(fit2)  # Estimated pA and f
#>        pA         f 
#> 0.5403226 0.0000000 
summary(fit1)
#> 
#> Call:
#> vglm(formula = y ~ 1, family = AA.Aa.aa, trace = TRUE)
#> 
#> Coefficients: 
#>             Estimate Std. Error z value Pr(>|z|)
#> (Intercept)   0.1616     0.1040   1.554     0.12
#> 
#> Name of linear predictor: logitlink(pA) 
#> 
#> Residual deviance: 0.1479 on 0 degrees of freedom
#> 
#> Log-likelihood: -5.384 on 0 degrees of freedom
#> 
#> Number of Fisher scoring iterations: 3 
#> 
#> No Hauck-Donner effect found in any of the estimates
#>