| Title: | Parental Allele Transmission Inference for Trio Genotype Data |
|---|---|
| Description: | Infers maternal and paternal transmitted and non-transmitted alleles from phased trio genotype data. The package supports SNP-level analyses of genetic nurture and transgenerational effects. It interoperates with Bioconductor VCF infrastructure through support for VariantAnnotation::VCF objects and returns R objects for downstream analysis. |
| Authors: | Jinyi Che [aut, cre] |
| Maintainer: | Jinyi Che <[email protected]> |
| License: | GPL-3 + file LICENSE |
| Version: | 1.1.0 |
| Built: | 2026-05-30 09:37:45 UTC |
| Source: | https://github.com/bioc/parati |
Given trio genotype data for a single family, infer maternal and paternal transmitted and non-transmitted alleles.
haplotype_infer(vcf_dt, hap_length = 5e+05)haplotype_infer(vcf_dt, hap_length = 5e+05)
vcf_dt |
A 'data.table' containing fixed VCF columns and genotype columns named 'M', 'P', and 'B'. |
hap_length |
Integer, haplotype window length. |
This function preserves the original PARATI inference semantics: - deterministic Mendelian inference for non-triple-heterozygote patterns - local haplotype matching for triple-heterozygote sites - low-similarity / ambiguous sites are left as missing (".|.")
A named list containing:
A 'data.table' with transmitted alleles.
A 'data.table' with non-transmitted alleles.
A 'data.table' summarizing inference statistics.
Main function to infer parental transmitted and non-transmitted alleles from phased trio genotype data. The function accepts either a VCF file path or a 'VCF' object from 'VariantAnnotation', and returns R objects without writing files by default.
parati_run(vcf, fam, chr = NULL, hap_length = 5e+05)parati_run(vcf, fam, chr = NULL, hap_length = 5e+05)
vcf |
Either a character path to a phased VCF/VCF.GZ file, or a 'VariantAnnotation::VCF' object. |
fam |
Either a character path to a family index table ('.xlsx') or a 'data.frame' / 'data.table' with columns 'FamilyIndex', 'IndividualID', and 'Role'. |
chr |
Optional chromosome identifier to subset variants. |
hap_length |
Integer, haplotype window length. |
A named list containing:
A 'data.table' with transmitted alleles.
A 'data.table' with non-transmitted alleles.
A 'data.table' summarizing inference statistics.
vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") fam_file <- system.file("extdata", "Toy_FamilyIndexTable.xlsx", package = "parati") res <- parati_run(vcf = vcf_file, fam = fam_file, chr = 1) names(res)vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") fam_file <- system.file("extdata", "Toy_FamilyIndexTable.xlsx", package = "parati") res <- parati_run(vcf = vcf_file, fam = fam_file, chr = 1) names(res)
Reads a VCF file and subsets variants to a given chromosome.
read_vcf_by_chr(vcf_file, chr)read_vcf_by_chr(vcf_file, chr)
vcf_file |
Character scalar, path to a VCF/VCF.GZ file. |
chr |
Character or integer chromosome identifier. |
A 'data.table' containing VCF rows for the selected chromosome.
vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_chr <- read_vcf_by_chr(vcf_file, chr = 1) dim(vcf_chr)vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_chr <- read_vcf_by_chr(vcf_file, chr = 1) dim(vcf_chr)
Converts a 'data.table' representing VCF rows into a 'vcfR' object.
vcf_dt_to_vcfR(df, meta = character())vcf_dt_to_vcfR(df, meta = character())
df |
A 'data.table' containing VCF rows. |
meta |
Character vector of VCF meta lines. |
A 'vcfR' object.
vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_dt <- read_vcf_by_chr(vcf_file, chr = 1) vcf_obj <- vcf_dt_to_vcfR(vcf_dt) class(vcf_obj) stopifnot(inherits(vcf_obj, "vcfR"))vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_dt <- read_vcf_by_chr(vcf_file, chr = 1) vcf_obj <- vcf_dt_to_vcfR(vcf_dt) class(vcf_obj) stopifnot(inherits(vcf_obj, "vcfR"))
Writes a 'data.table' representing VCF rows to a VCF file.
write_vcf_dt(df, file)write_vcf_dt(df, file)
df |
A 'data.table' containing VCF rows. |
file |
Character scalar, output file path. |
'NULL', invisibly. The VCF file is written to 'file'.
vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_dt <- read_vcf_by_chr(vcf_file, chr = 1) outfile <- tempfile(fileext = ".vcf") write_vcf_dt(vcf_dt, outfile) file.exists(outfile)vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_dt <- read_vcf_by_chr(vcf_file, chr = 1) outfile <- tempfile(fileext = ".vcf") write_vcf_dt(vcf_dt, outfile) file.exists(outfile)
Writes a 'vcfR' object to a VCF file.
write_vcf_obj(vcf_obj, file)write_vcf_obj(vcf_obj, file)
vcf_obj |
A 'vcfR' object. |
file |
Character scalar, output file path. |
'NULL', invisibly. The VCF file is written to 'file'.
vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_dt <- read_vcf_by_chr(vcf_file, chr = 1) vcf_obj <- vcf_dt_to_vcfR(vcf_dt) stopifnot(inherits(vcf_obj, "vcfR")) outfile <- tempfile(fileext = ".vcf") write_vcf_obj(vcf_obj, outfile) stopifnot(file.exists(outfile))vcf_file <- system.file("extdata", "Toy_TrioGenotype.vcf.gz", package = "parati") vcf_dt <- read_vcf_by_chr(vcf_file, chr = 1) vcf_obj <- vcf_dt_to_vcfR(vcf_dt) stopifnot(inherits(vcf_obj, "vcfR")) outfile <- tempfile(fileext = ".vcf") write_vcf_obj(vcf_obj, outfile) stopifnot(file.exists(outfile))