R/RcppExports.R

In eriqande/fullsniplings: Full sibling reconstruction from SNP data

# This file was generated by Rcpp::compileAttributes
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393

#' Compute the likelihood of the genotypes of a marriage given the genotype of \emph{just one} of its offspring
#' 
#' It is advantageous to precompute and store some quantities for pedigree analysis with SNPs.
#' One of those quantities is the likelihood of the genotypes of the two parents that produced a child,
#' conditional only on that one child's genotype. This is convenient because the likelihood of the two parents' 
#' genotypes given \emph{all} of their offspring will be the product of the single-offspring
#' likelihoods.  When genotyping error is present, that has to be thrown into the mix as well.
#' However, here we assume that the individual offspring genotype likelihoods have already been
#' computed given their observed genotypes.
#' @param offspring_likelihoods This is a 3 x L x N array of likelihoods of the individual's
#' true genotypes given their observed genotypes.  This sort of array is returned by the function
#' \code{\link{get_indiv_geno_lik}}.  
#' @param transmission_probs a 3 x 3 x 3 matrix of probs of kid genotypes given parent genotypes.
#' Such an array is returned by \code{\link{trans_probs}}.  See the documentation thereof for a
#' description.
#' @return This returns a 3 x 3 x L x N array of the parent genotypes (9 states) at each of the L loci
#' for each of the N offspring.  Although these are ostensible likelihoods of "marriages" (think of them
#' as marriage nodes), it will be beneficial to regard these quantities of properties of the individual
#' offspring moreso than as properties of the marriage nodes.
#' @export
per_kid_marriage_likelihoods <- function(offspring_likelihoods, transmission_probs) {
    .Call('fullsniplings_per_kid_marriage_likelihoods', PACKAGE = 'fullsniplings', offspring_likelihoods, transmission_probs)
}

#' compute the genotype probs of a pair of full siblings with genotyping error
#' 
#' This could be done entirely in R, but it is hard to think it through and it seemed
#' it would be hard to maintain.  Super easy in C.
#' 
#' @param  L the number of loci
#' @param G the number of genotypic states
#' @param PP parent pair probs as returned by unrelated_pair_gfreqs() for example
#' @param TP transmision probs as returned by trans_probs()
#' @param GE prob of observed genotypes given true genotypes, as returned by lik_array_from_simple_geno_err() for example
#' @export
C_full_sibling_pair_gfreqs <- function(L, G, PP, TP, GE) {
    .Call('fullsniplings_C_full_sibling_pair_gfreqs', PACKAGE = 'fullsniplings', L, G, PP, TP, GE)
}

#' compute the marriage node likelihoods given offspring specified in a full-sibling list
#' 
#' @param S a list of vectors that give the indices (base 0) of the individuals in the 
#' full sibling groups.  
#' @param PK per-kid marriage likelihoods.  This must be of class \code{\link{marriage_geno_lik_array}},
#' which is just a matrix underneath with G x G x L rows and N columns.
#' @export
multi_kid_marriage_likelihoods <- function(S, PK) {
    .Call('fullsniplings_multi_kid_marriage_likelihoods', PACKAGE = 'fullsniplings', S, PK)
}

#' Update one row of a marriage likelihood matrix IN PLACE!
#' 
#' The intended use of this is to update the marriage likelihoods after an individual has
#' been moved from one sibship to another.  So, for example, if you moved an individual 
#' from sibship 40 to sibship 291 (as subscripted in R) then you would want to run this
#' with \code{bz_idx} equal to \code{c(39, 290)}.  
#' @param S a list of vectors that give the indices (base 0) of the individuals in the 
#' full sibling groups.  
#' @param PK per-kid marriage likelihoods.  This must be of class \code{\link{marriage_geno_lik_array}},
#' which is just a matrix underneath with G x G x L rows and N columns.
#' @param ML the marriage likelihoods matrix to be modified
#' @param bz_idx An integer vector holding the BASE-0 indices of the components of S that 
#' will be accessed and used to update ML.
#' 
#' @return This doesn't return anything.  It modifies ML in place via call be reference.  Our
#' goal here is to make updates without copying a lot of memory.
#' @export
update_marriage_likelihoods_in_place <- function(S, PK, ML, bz_idx) {
    invisible(.Call('fullsniplings_update_marriage_likelihoods_in_place', PACKAGE = 'fullsniplings', S, PK, ML, bz_idx))
}

#' Calculate/update one row of a marriage Posterior matrix IN PLACE!
#' 
#' The intended use of this is to update the marriage posteriors after an individual has
#' been moved from one sibship to another.  So, for example, if you moved an individual 
#' from sibship 40 to sibship 291 (as subscripted in R) then you would want to run this
#' with \code{bz_idx} equal to the column of sibship 39, and then run it again with it equal 
#' to the column of  290, after you had updated the LMMFS
#' for each of those.
#' @param ML the marriage likelihoods matrix
#' @param MP the marriage posteriors matrix that is to be updated
#' @param Pri  the joint prior probbabilities of the parent pair
#' @param NGS Number of genotypic states.  For pairs of parents, for example, this will be 9. 
#' @param bz_idx An integer vector holding the BASE-0 index of the columns of MP that needs updating.  Note that this
#' is not the index of the sibhip, but rather the column that sibship occupies in MP.  That will be something like 
#' FSL[Indiv]$LMMI_Idx that 
#' 
#' @return This returns the value of the index of the last marriage that it updated.  It does this so I can check some things.  It modifies MP in place via call by reference.  Our
#' goal here is to make updates without copying a lot of memory.
#' @export
update_marriage_posteriors_in_place <- function(ML, MP, Pri, NGS, bz_idx) {
    .Call('fullsniplings_update_marriage_posteriors_in_place', PACKAGE = 'fullsniplings', ML, MP, Pri, NGS, bz_idx)
}

#' update marriage node kid prongs IN PLACE
#' 
#' This uses the information in the PMMFS matrix (the 
#' "Posterior Matrix of Marriages given Full Siblings") to calculate
#' the posterior predictive distribution for the next full sibling
#' in the full sibling group which is referenced through S and bz_idx
#'  
#' @param S a list of lists each with two components.  The first is LMMI_Idx which is the base-0
#' index of the Marriage that the component is referring to, and the second is Indivs which give
#' the indices (base 0) of the individuals in the 
#' full sibling groups.  
#' @param MP the marriage posteriors matrix to be used in the udpating. Should be a matrix
#' that has NGS_P * L rows, where L is the number of loci.
#' @param KP Kid prongs.  This is what gets updated in place. Should be a matrix that has
#' NGS_K * L rows, where L is the number of loci.
#' @param NGS_P Number of genotypic states in a parent pair.  This will typically be 9.
#' @param NGS_K Number of genotypic states in a kid.  This will typically be 3
#' @param TP transmision probs as returned by trans_probs()
#' @param bz_idx An integer vector holding the BASE-0 indices of the components of S that 
#' will be accessed and used to update ML.
#' 
#' @return This doesn't return anything.  It modifies Prongs in place via call by reference.  Our
#' goal here is to make updates without copying a lot of memory.
update_marriage_node_kid_prongs_in_place <- function(S, MP, KP, NGS_P, NGS_K, TP, bz_idx) {
    invisible(.Call('fullsniplings_update_marriage_node_kid_prongs_in_place', PACKAGE = 'fullsniplings', S, MP, KP, NGS_P, NGS_K, TP, bz_idx))
}

#' Computes the posterior predictive likelihoods that an individual with likelihood vector IndGenoLik belongs to each sibling group in AFS
#' 
#' Note that we have to pass in the whole full sibling list (FSL) so that we can get the LMMI_Idx for each sibling group.
#' @param IndGenoLik Likelihoods of the focal individual.  There will be 3*L elements in that.  Note that \code{AFS} is the acceptable
#' sibling list for this focal individual.
#' @export
kid_prongs_times_ind_likelihoods <- function(FSL, IndGenoLik, KidProngs, AFS) {
    .Call('fullsniplings_kid_prongs_times_ind_likelihoods', PACKAGE = 'fullsniplings', FSL, IndGenoLik, KidProngs, AFS)
}

#' return a vector of indices of sibships that include at least one hi-sibship-lod individual
#' 
#' Note that it is currently up to the user to ensure that no element of AFS exceeds the length-1
#' of IFS
#' 
#' @param IFS the individual full siblings group vector. A vector such that element i contains the 
#' base-0 index of the sibship to which the individual with index i belongs.
#' @param AFS the vector of individual indexes (base 0) that have hi sibship lod with the 
#' focal individual.  i.e. AFS = the Acceptable Full Siblings
#' @export
#' @examples
#' possible_sibgroups(rep(16,50), 0:49)
#' possible_sibgroups(rep(0:9, each=3), sample(0:29, 10))
possible_sibgroups <- function(IFS, AFS) {
    .Call('fullsniplings_possible_sibgroups', PACKAGE = 'fullsniplings', IFS, AFS)
}

#' compute the pseudo-prior for individual and individual pulled out of sibgroup IndG
#' 
#' @param FSL the full sibling list
#' @param IndG The full sibling group to which the focal individual belongs
#' @param AFS acceptable full sibling groups for the focal individual.
#' @return This returns a list of two components. $solo is the prob that the individual will form a singleton and 
#' $afs is the prior that the individual will join any of the existing sibgroups in afs, scaled to sum to one.
#' @export
pseudo_prior <- function(FSL, IndG, AFS) {
    .Call('fullsniplings_pseudo_prior', PACKAGE = 'fullsniplings', FSL, IndG, AFS)
}

#' compute the posterior of Inds genotype given the genotype freqs in the pop
#' 
#' @param Gfreqs  array of genotype frequencies
#' @param Liks  array of likleihoods (3 * L in length)
#' @export
geno_post_c <- function(Gfreqs, Liks) {
    .Call('fullsniplings_geno_post_c', PACKAGE = 'fullsniplings', Gfreqs, Liks)
}

#' does a gibbs update of the full sibling group of individual Ind
#' @param FSL list of integer vectors. Each component is a list of two componontes: \code{LMMI_Idx} = the index of the sibgroup
#' in LMMI (see below).  It also is the column that Ind's current full sib group occuppies in LMMFS
#' @param IFS individual full siblings group vector. A vector such that element i contains the 
#' base-0 index of the sibship to which the individual with index i belongs.
#' and PMMFS, and KidProngs \code{Indivs} = the the base-0 indices of individuals in the full sibship.
#' @param LMMI Likelihood matrix of marriages given individuals.
#' @param LMMFS Likelihood matrix of marriages given full sibships.
#' @param PMMFS Posterior matrix of marriages given full sibships.
#' @param KidProngs Posterior predictives for the next individuals to be sampled from a full sibship.
#' @param Pile  An integer vector that we will use as a stack to hold the indices of empty FSL list elements
#' that we can put newly formed sibships into
#' @param MatPile An integer vector that is parallel to Pile that holds the corresponding column of LMMFS that 
#' go along with the elements of Pile. 
#' can be renewed.
#' @param AFSL The acceptable full siblings list.  Component i is a vector of the base-0 indices of the individuals
#' that indiv i has high full-sibling LOD with.
#' @param Gfreqs  The genotype frequencies.
#' @param UPG Unrelated Pair Genotype Frequencies.
#' @param TP Transmission probabilities
#' @param IndLiks  Individual likelihoods
#' @param Ind the index of the individual to be updated.  
#' 
#' 
#' @export
gibbs_update_one_indiv_in_place <- function(FSL, IFS, LMMI, LMMFS, PMMFS, KidProngs, Pile, MatPile, AFSL, Gfreqs, UPG, TP, IndLiks, Ind) {
    .Call('fullsniplings_gibbs_update_one_indiv_in_place', PACKAGE = 'fullsniplings', FSL, IFS, LMMI, LMMFS, PMMFS, KidProngs, Pile, MatPile, AFSL, Gfreqs, UPG, TP, IndLiks, Ind)
}

#' return the high-sib-pair-logl indivs for each indiv
#' 
#' I wrote this to see how much faster I could implement the functionality of find_high_logl_sib_pairs
#' using Rcpp.  find_high_logl_sib_pairs is painfully slow, and I think this should be much faster
#' @param FSP a 9 x L matrix (3 x 3 x L) of expected genotype frequencies of a full sibling pair,
#' but we pass it in is a vector.
#' @param UPF a 9 x L matrix of expected genotype frequencies of an unrelated pair, which we also
#' pass in as a vector.
#' @param G an L x N matrix of 0, 1, or 2 or NA giving the genotypes of the individuals 
#' @param loglV the cutoff point above which you will accept the pairs.
#' @export
high_logl_pairs <- function(FSP, UPF, G, loglV) {
    .Call('fullsniplings_high_logl_pairs', PACKAGE = 'fullsniplings', FSP, UPF, G, loglV)
}