Package {BIGpopA}

Compute Allele Frequencies for Populations

Description

Computes allele frequencies for specified populations given SNP array data.

Usage

allele_freq_poly(geno, populations, ploidy = 2)

Arguments

Details

When geno is a .ped file, the returned matrix carries a counted_allele attribute. solve_composition_poly() uses it to code a validation .ped file with the same counted allele at each marker.

Value

A matrix of allele frequencies with SNPs in rows and populations in columns.

References

Funkhouser SA, Bates RO, Ernst CW, Newcom D, Steibel JP. Estimation of genome-wide and locus-specific breed composition in pigs. Transl Anim Sci. 2017 Feb 1;1(1):36-44.

Examples

geno_matrix <- matrix(
  c(4, 1, 4, 0,
    2, 2, 1, 3,
    0, 4, 0, 4,
    3, 3, 2, 2,
    1, 4, 2, 3),
  nrow = 4, ncol = 5, byrow = FALSE,
  dimnames = list(paste0("Ind", 1:4), paste0("S", 1:5))
)

pop_list <- list(
  PopA = c("Ind1", "Ind2"),
  PopB = c("Ind3", "Ind4")
)

allele_freqs <- allele_freq_poly(geno = geno_matrix,
                                 populations = pop_list,
                                 ploidy = 4)
print(allele_freqs)

Check and Correct Common Pedigree Errors

Description

Reads a 3-column pedigree file (id, male_parent, female_parent) and performs quality checks, optionally correcting detected errors. Exact duplicates and missing parents are always corrected. Conflicting trios and inconsistent sex roles are corrected when their respective arguments are TRUE. Cycles are reported only and must be resolved manually.

Usage

check_ped(
  ped.file,
  seed = NULL,
  verbose = TRUE,
  correct_conflicting_trios = TRUE,
  correct_inconsistent_sex_roles = TRUE
)

Arguments

Value

An invisible named list of data frames:

exact_duplicates

Exact duplicate rows found in the input.

conflicting_trios

IDs with conflicting male_parent or female_parent assignments.

inconsistent_sex_roles

Rows where a conflicting ID appears as male_parent or female_parent.

missing_parents

Parent IDs absent from id, added as founders.

dependencies

Cycles detected in the pedigree. Must be resolved manually.

corrected_pedigree

Corrected pedigree table.

Author(s)

Josue Chinchilla-Vargas

Examples

# Self-contained example using a data.frame
ped_df <- data.frame(
  id            = c("A", "B", "C", "C", "D"),
  male_parent   = c("0", "0", "A", "A", "B"),
  female_parent = c("0", "0", "B", "B", "C"),
  stringsAsFactors = FALSE
)
ped_errors <- check_ped(ped.file = ped_df, seed = 101919, verbose = FALSE)
names(ped_errors)
head(ped_errors$corrected_pedigree)


library(data.table)
ped_dt <- data.table(id = c("A", "B", "C"),
                     male_parent   = c("0", "0", "A"),
                     female_parent = c("0", "0", "B"))
ped_errors <- check_ped(ped.file = ped_dt, verbose = FALSE)

Find Parentage Assignments for Progeny

Description

Assigns the most likely parent(s) to each progeny from SNP genotype data using Mendelian error rates or homozygous mismatch rates. Parents or progeny absent from the genotype file are removed with a warning.

Usage

find_parentage(
  genotypes_file,
  parents_file,
  progeny_file,
  method = "best_pair",
  min_markers = 10,
  error_threshold = 5,
  show_ties = TRUE,
  allow_parent_selfing = FALSE,
  exclude_self_match = TRUE,
  verbose = TRUE,
  plot_results = TRUE,
  ploidy = 2
)

Arguments

Value

A named list (returned invisibly) with elements:

pass

Progeny with a confident parentage assignment.

high_error

Progeny whose best assignment exceeds the error threshold.

low_markers

Progeny with insufficient markers for a valid assignment.

full_results

Complete data.table with all progeny and all output columns.

plot

ggplot object if plot_results = TRUE, otherwise NULL.

Author(s)

Josue Chinchilla-Vargas

Examples

geno_df <- data.frame(
  id  = c("P1", "P2", "P3", "Off1", "Off2"),
  S1  = c(0L, 2L, 0L, 1L, 0L),
  S2  = c(2L, 0L, 2L, 1L, 2L),
  S3  = c(0L, 2L, 0L, 1L, 0L),
  S4  = c(2L, 0L, 2L, 1L, 2L),
  S5  = c(0L, 2L, 0L, 1L, 0L),
  S6  = c(2L, 0L, 2L, 1L, 2L),
  S7  = c(0L, 2L, 0L, 1L, 0L),
  S8  = c(2L, 0L, 2L, 1L, 2L),
  S9  = c(0L, 2L, 0L, 1L, 0L),
  S10 = c(2L, 0L, 2L, 1L, 2L)
)

parents_df <- data.frame(
  id  = c("P1", "P2", "P3"),
  sex = c("M",  "F",  "F"),
  stringsAsFactors = FALSE
)

progeny_df <- data.frame(
  id = c("Off1", "Off2"),
  stringsAsFactors = FALSE
)

results <- find_parentage(
  genotypes_file = geno_df,
  parents_file   = parents_df,
  progeny_file   = progeny_df,
  method         = "best_pair",
  verbose        = FALSE,
  plot_results   = FALSE
)
print(results$full_results)

Convert a PLINK .ped File to a BIGpopA Genotype Matrix

Description

Reads a PLINK 1 text pedigree file (.ped) and returns allele dosages in the layout used by the rest of BIGpopA: individuals in rows, markers in columns, coded as the number of copies of the counted allele (0, 1, 2) with NA for missing calls. PLINK .ped files are diploid only.

Usage

ped_to_dosage(
  ped,
  map = "auto",
  format = c("data.frame", "matrix"),
  counted_allele = NULL,
  missing_codes = c("0", "N", "-", "."),
  verbose = TRUE
)

Arguments

Details

The .ped format stores allele letters rather than a reference/alternate pair, so the counted allele (allele B) is chosen per marker. By default it is the alphabetically / numerically last allele observed at the marker, so 1/2 coded files count allele 2 and an A/G marker counts G. To code a second file the same way as a first one (for example validation genotypes against a reference panel), pass the first file's counted_allele attribute through the counted_allele argument. Otherwise a marker that is monomorphic in one file can be coded in the opposite direction.

Only the individual ID (column 2) is used from the six leading pedigree columns. A .map file is optional and only supplies marker names; it is not used for allele coding. Without it, markers are named SNP1, SNP2, ... in file order.

Value

A data.table with an id column followed by marker columns (format = "data.frame"), or an integer matrix of individuals by markers (format = "matrix"). The attribute "ploidy" is 2 and the attribute "counted_allele" holds the allele counted at each returned marker. Markers with more than two alleles are dropped with a warning.

Author(s)

Josue Chinchilla-Vargas

See Also

vcf_to_dosage(), find_parentage(), allele_freq_poly()

Examples

ped_path <- tempfile(fileext = ".ped")
writeLines(c(
  "F1 P1 0 0 1 -9 A A G G A G",
  "F1 P2 0 0 2 -9 A G G G G G",
  "F1 O1 P1 P2 1 -9 A G G G A G"
), ped_path)

genotypes <- ped_to_dosage(ped_path, map = NULL, verbose = FALSE)
print(genotypes)
attr(genotypes, "counted_allele")

unlink(ped_path)

Compute Genome-Wide Breed Composition

Description

Computes genome-wide breed/ancestry composition using quadratic programming on a batch of animals.

Usage

solve_composition_poly(Y, X, ploidy = 2)

Arguments

Value

A matrix with one row per animal, one column per reference population (estimated proportions, summing to 1), and an R2 column.

References

Funkhouser SA, Bates RO, Ernst CW, Newcom D, Steibel JP. Estimation of genome-wide and locus-specific breed composition in pigs. Transl Anim Sci. 2017 Feb 1;1(1):36-44.

Examples

allele_freqs_matrix <- matrix(
  c(0.625, 0.500,
    0.500, 0.500,
    0.500, 0.500,
    0.750, 0.500,
    0.625, 0.625),
  nrow = 5, ncol = 2, byrow = TRUE,
  dimnames = list(paste0("SNP", 1:5), c("VarA", "VarB"))
)

val_geno_matrix <- matrix(
  c(2, 1, 2, 3, 4,
    3, 4, 2, 3, 0),
  nrow = 2, ncol = 5, byrow = TRUE,
  dimnames = list(paste0("Test", 1:2), paste0("SNP", 1:5))
)

composition <- solve_composition_poly(Y = val_geno_matrix,
                                      X = allele_freqs_matrix,
                                      ploidy = 4)
print(composition)

Validate Pedigree Trios Using Mendelian Error Analysis

Description

Validates parent-offspring trios against SNP genotype data using Mendelian error rates. Identifies incorrect parentage assignments, suggests best-matching replacements, and outputs a corrected pedigree. Founder trios (both parents coded as 0) are preserved unchanged if a founders file is supplied. Trios absent from the genotype file are retained as no_genotype_data.

Usage

validate_pedigree(
  pedigree_file,
  genotypes_file,
  founders_file = NULL,
  trio_error_threshold = 5,
  min_markers = 10,
  single_parent_error_threshold = 2,
  verbose = TRUE,
  plot_results = TRUE,
  ploidy = 2
)

Arguments

Value

An invisible named list with the following elements:

pass

Trios that passed the Mendelian error threshold.

fail

Trios that failed the Mendelian error threshold.

low_markers

Trios with insufficient markers for evaluation.

no_genotype_data

Trios absent from the genotype file.

founders

Trios identified as founders.

missing_parents

Trios with one or both parents coded as 0 (non-founders).

full_results

Complete data.table with all trios and all output columns.

corrected_pedigree

Pedigree table after applying recommended corrections.

plot

ggplot object if plot_results = TRUE, otherwise NULL.

Author(s)

Josue Chinchilla-Vargas

Examples

geno_df <- data.frame(
  id  = c("P1", "P2", "P3", "Off1", "Off2"),
  S1  = c(0L, 2L, 0L, 1L, 0L),
  S2  = c(2L, 0L, 2L, 1L, 2L),
  S3  = c(0L, 2L, 0L, 1L, 0L),
  S4  = c(2L, 0L, 2L, 1L, 2L),
  S5  = c(0L, 2L, 0L, 1L, 0L),
  S6  = c(2L, 0L, 2L, 1L, 2L),
  S7  = c(0L, 2L, 0L, 1L, 0L),
  S8  = c(2L, 0L, 2L, 1L, 2L),
  S9  = c(0L, 2L, 0L, 1L, 0L),
  S10 = c(2L, 0L, 2L, 1L, 2L)
)

ped_df <- data.frame(
  id            = c("Off1", "Off2"),
  male_parent   = c("P1",   "P1"),
  female_parent = c("P2",   "P3"),
  stringsAsFactors = FALSE
)

results <- validate_pedigree(
  pedigree_file  = ped_df,
  genotypes_file = geno_df,
  verbose        = FALSE,
  plot_results   = FALSE
)
print(results$full_results)

Convert a VCF File to a BIGpopA Genotype Matrix

Description

Reads a variant call format (VCF) file, or an existing vcfR object, and returns allele dosages in the layout used by the rest of BIGpopA: individuals in rows, markers in columns, coded as the dosage of the alternate (B) allele (0, 1, ..., ploidy) with NA for missing calls.

Usage

vcf_to_dosage(
  vcf,
  ploidy = NULL,
  format = c("data.frame", "matrix"),
  dosage_field = NULL,
  ref_dosage = FALSE,
  biallelic_only = TRUE,
  marker_names = c("auto", "id", "position"),
  min_marker_call_rate = 0,
  min_sample_call_rate = 0,
  verbose = TRUE,
  ...
)

Arguments

Details

The function is ploidy agnostic. Dosage is obtained by counting the non-reference alleles in each GT call, so a diploid 0/1 returns 1 and a tetraploid 0/1/1/1 returns 3, with no assumption about the ploidy level. When ploidy is not supplied it is inferred from the number of allele slots in the calls themselves (the most common slot count across the file). Calls whose slot count disagrees with the file ploidy, for instance a diploid call in an otherwise tetraploid file, are set to missing and reported, since the downstream functions require a single ploidy level.

Phased (|) and unphased (/) separators are both accepted. A call is treated as missing when any of its allele slots is ..

vcfR is used to read the file and is listed under Suggests, so it must be installed before calling this function (install.packages("vcfR")).

Value

A data.table with an id column followed by marker columns (format = "data.frame"), or a numeric matrix of individuals by markers (format = "matrix"). The attribute "ploidy" records the ploidy used.

Author(s)

Josue Chinchilla-Vargas

See Also

find_parentage(), validate_pedigree(), allele_freq_poly()

Examples

if (requireNamespace("vcfR", quietly = TRUE)) {
  # Write a small tetraploid VCF (3 markers x 4 samples) to a temporary file
  vcf_path <- tempfile(fileext = ".vcf")
  writeLines(c(
    "##fileformat=VCFv4.3",
    "##FORMAT=<ID=GT,Number=1,Type=String,Description=\"Genotype\">",
    paste("#CHROM", "POS", "ID", "REF", "ALT", "QUAL", "FILTER", "INFO",
          "FORMAT", "S1", "S2", "S3", "S4", sep = "\t"),
    paste("chr1", "100", "snp1", "A", "T", ".", ".", ".", "GT",
          "0/0/0/0", "0/0/0/1", "0/1/1/1", "1/1/1/1", sep = "\t"),
    paste("chr1", "200", "snp2", "C", "G", ".", ".", ".", "GT",
          "0/0/1/1", "./././.", "0/0/0/0", "0/1/1/1", sep = "\t"),
    paste("chr1", "300", "snp3", "G", "A", ".", ".", ".", "GT",
          "1/1/1/1", "1/1/1/1", "0/0/1/1", "0/0/0/1", sep = "\t")
  ), vcf_path)

  # tetraploid VCF -> id + marker columns for find_parentage()
  genotypes <- vcf_to_dosage(vcf_path, ploidy = 4)
  print(genotypes)

  # same data as an individuals x markers matrix for allele_freq_poly()
  geno_mat <- vcf_to_dosage(vcf_path, ploidy = 4, format = "matrix")

  unlink(vcf_path)
}