fEstimate: Estimation of inbreeding coefficients

View source: R/fEstimate.R

fEstimateR Documentation

Estimation of inbreeding coefficients

Description

Estimate the inbreeding coefficient f of pedigree members from their genotypes. The default method uses maximum likelihood, but moment estimators are also available.

Usage

fEstimate(
  x,
  ids = typedMembers(x),
  method = c("mle", "simple", "ritland"),
  markers = NULL,
  verbose = FALSE
)

Arguments

x

A ped object or a list of such.

ids

A vector with ID labels, by default, all genotyped members of x.

method

A character, either "mle" (maximum-likelihood estimate; default), "simple" (simple method-of-moments estimator), or "ritland" (see Details).

markers

A vector with names or indices of markers attached to x, indicating which markers to include. By default, all markers are included. Missing genotypes are removed before estimation.

verbose

A logical.

Details

The ML estimator assumes independent markers and known allele frequencies. For a single marker with allele frequencies p_a, the genotype likelihood P(G \mid f) is:

P(G = a/a \mid f) = f p_a + (1-f)p_a^2,

P(G = a/b \mid f) = 2(1-f)p_a p_b,\quad a \neq b.

For multiple markers, the full likelihood is the product over markers. The function uses optimise() to maximise the log-likelihood over [0,1].

The "simple" method is a method-of-moments estimator based on observed and expected homozygosity. Let L be the number of markers, N_H the number of observed homozygous genotypes, and H^* = \sum_1^L \sum_a p_a^2 the expected number under Hardy–Weinberg equilibrium. The estimator is then:

\hat{f} = \dfrac{N_H - H^*}{L - H^*}.

The "ritland" method implements the estimator of Ritland (1996). For marker j, let n_j be the number of alleles, and define T_j to be 1/p_a if the genotype is a/a, and 0 otherwise. The estimator is then:

\hat{f} = \dfrac{\sum_{j=1}^L (T_j - 1)}{\sum_{j=1}^L (n_j - 1)} = \dfrac{\sum_{\text{hom}} 1/p_a - L}{A - L},

where A = \sum n_j is the total number of alleles across markers. Note that both this and the "simple" estimator may produce estimates outside the interval [0,1].

The MLE approach is described explicitly in Section 3.2 of Thompson (1986). Ritland (1996) introduces the Ritland estimator, and also mentions both the simple estimator and the MLE. Note that Ritland's formulation of the MLE is slightly different (but equivalent) to Thompson's.

Value

A named numeric vector.

References

  • Thompson, E. A. (1986). Pedigree Analysis in Human Genetics. Johns Hopkins University Press.

  • Ritland, K. (1996). Estimators for pairwise relatedness and individual inbreeding coefficients. Genetical Research.

See Also

ibdEstimate()

Examples


# Simulate genotype for 100 SNPs, for a child of full sibs
x = nuclearPed(2, sex = 1:2) |> addSon(3:4) |>
  markerSim(N = 100, ids = 5, alleles = 1:2, afreq = c(0.3, 0.7), seed = 123)
x

# Estimate inbreeding coefficient (pedigree: 0.25)
fEstimate(x) # MLE
fEstimate(x, method = "simple")
fEstimate(x, method = "ritland")


#-----------------------------
# Compare different estimators
#-----------------------------

# Simulate 200 DNA profiles (35 STR markers) for a child of half-sibs
x = halfSibPed(sex2 = 2) |> addSon(4:5) |>
  profileSim(N = 200, ids = 6, markers = NorwegianFrequencies, seed = 123)

# Estimates using all three methods
fhat = data.frame(
  mle = sapply(x, fEstimate, method = "mle"),
  sim = sapply(x, fEstimate, method = "simple"),
  rit = sapply(x, fEstimate, method = "ritland")
)

# Mean estimates (expected: 0.125)
apply(fhat, 2, mean)

# RMSE
apply(fhat, 2, function(v) sqrt(mean((v - 0.125)^2)))


forrel documentation built on June 28, 2026, 5:07 p.m.