From 2782d6dcc9ee4be9855b5e468ce789425b81d49a Mon Sep 17 00:00:00 2001 From: gilesjohnr Date: Mon, 17 Dec 2018 13:00:34 -0500 Subject: [PATCH] fixing half-closed notation and examples --- DESCRIPTION | 8 +- NAMESPACE | 5 + R/examples/get_cross_k.R | 3 +- R/examples/get_tau.R | 24 +- R/spatialfuncs.r | 86 +- R/transdistfuncs.r | 49 +- R/wrapperfuncs.r | 52 +- man/est.transdist.Rd | 26 +- man/est.transdist.bootstrap.ci.Rd | 30 +- man/est.transdist.temporal.Rd | 43 +- man/est.transdist.temporal.bootstrap.ci.Rd | 55 +- man/est.transdist.theta.weights.Rd | 41 +- man/est.wt.matrix.Rd | 8 +- man/est.wt.matrix.weights.Rd | 10 +- man/get.cross.K.Rd | 7 +- man/get.cross.PCF.Rd | 4 +- man/get.pi.Rd | 26 +- man/get.pi.bootstrap.Rd | 7 +- man/get.pi.ci.Rd | 11 +- man/get.pi.permute.Rd | 10 +- man/get.pi.typed.Rd | 13 +- man/get.pi.typed.bootstrap.Rd | 7 +- man/get.pi.typed.permute.Rd | 7 +- man/get.tau.Rd | 12 +- man/get.tau.bootstrap.Rd | 11 +- man/get.tau.ci.Rd | 11 +- man/get.tau.permute.Rd | 11 +- man/get.tau.typed.Rd | 20 +- man/get.tau.typed.bootstrap.Rd | 10 +- man/get.tau.typed.permute.Rd | 7 +- man/get.theta.Rd | 18 +- man/get.theta.bootstrap.Rd | 7 +- man/get.theta.ci.Rd | 7 +- man/get.theta.permute.Rd | 7 +- man/get.theta.typed.Rd | 17 +- man/get.theta.typed.bootstrap.Rd | 7 +- man/get.theta.typed.permute.Rd | 7 +- man/get.transdist.theta.Rd | 17 +- man/sim.epidemic.Rd | 34 +- src/spatialfuncs.c | 1040 ++++++++++---------- 40 files changed, 724 insertions(+), 1051 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 567b25d..781b622 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,14 +1,14 @@ Package: IDSpatialStats -Version: 0.3.5 -Date: 2018-06-11 +Version: 0.3.7 +Date: 2018-12-17 Title: Estimate Global Clustering in Infectious Disease -Author: Justin Lessler , Henrik Salje , John Giles +Author: Justin Lessler , Henrik Salje , John Giles Maintainer: Justin Lessler License: GPL (>=2) Description: Implements various novel and standard clustering statistics and other analyses useful for understanding the spread of infectious disease. -RoxygenNote: 6.0.1 +RoxygenNote: 6.1.1 Encoding: UTF-8 Imports: igraph, spatstat Depends: doParallel, foreach, parallel, R (>= 2.10) diff --git a/NAMESPACE b/NAMESPACE index 447e0bd..c8558ff 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -8,7 +8,10 @@ importFrom("igraph", "shortest.paths") # for get.transdist.theta importFrom("igraph", "E") # for get.transdist.theta importFrom("spatstat", "bounding.box.xy") # for est.transdist importFrom("spatstat", "as.ppp") # for est.transdist +importFrom("spatstat", "ppp") # for est.transdist importFrom("spatstat", "crossdist") # for est.transdist +importFrom("spatstat", "Kcross") # for get.cross.K +importFrom("spatstat", "pcf") # for get.cross.PCF export(get.pi) export(get.pi.ci) export(get.pi.bootstrap) @@ -40,5 +43,7 @@ export(est.transdist) export(est.transdist.bootstrap.ci) export(est.transdist.temporal) export(est.transdist.temporal.bootstrap.ci) +export(get.cross.K) +export(get.cross.PCF) useDynLib(IDSpatialStats, .registration=TRUE, .fixes="C_") diff --git a/R/examples/get_cross_k.R b/R/examples/get_cross_k.R index adab143..1480590 100644 --- a/R/examples/get_cross_k.R +++ b/R/examples/get_cross_k.R @@ -1,7 +1,6 @@ data(DengueSimulationR01) k <- get.cross.K(epi.data=DengueSimR01, type=5, hom=2, het=NULL, r=NULL, correction='border') -k.theo <- pi*(k$r)^2 -plot(k.theo, type='l', col='red', lty=2, xlab='r', ylab='cross K function') +plot(k[,2], type='l', col='red', lty=2, xlab='r', ylab='cross K function') lines(k$border) \ No newline at end of file diff --git a/R/examples/get_tau.R b/R/examples/get_tau.R index 9302b9e..934d735 100644 --- a/R/examples/get_tau.R +++ b/R/examples/get_tau.R @@ -7,21 +7,21 @@ r.min<-seq(0,980,20) r.mid<-(r.max+r.min)/2 sero.type.func<-function(a,b,tlimit=20){ - if(a[5]==b[5]&(abs(a[3]-b[3])<=tlimit)){rc=1} - else{rc=2} - return(rc) + if(a[5]==b[5]&(abs(a[3]-b[3])<=tlimit)){rc=1} + else{rc=2} + return(rc) } geno.type.func<-function(a,b,tlimit=20){ - if(a[4]==b[4]&(abs(a[3]-b[3])<=tlimit)){rc=1} - else{rc=2} - return(rc) + if(a[4]==b[4]&(abs(a[3]-b[3])<=tlimit)){rc=1} + else{rc=2} + return(rc) } sero.type.rep.func<-function(a,b,tlimit=20){ - if(a[5]==1&b[5]==1&(abs(a[3]-b[3])<=tlimit)){rc=1} - else{if(a[5]==1&b[5]==-999){rc=2}else{rc=3}} - return(rc) + if(a[5]==1&b[5]==1&(abs(a[3]-b[3])<=tlimit)){rc=1} + else{if(a[5]==1&b[5]==-999){rc=2}else{rc=3}} + return(rc) } sero.tau.R01<-get.tau(DengueSimR01,sero.type.func,r=r.max,r.low=r.min,comparison.type="independent") @@ -31,7 +31,7 @@ sero.tau.R02<-get.tau(DengueSimR02,sero.type.func,r=r.max,r.low=r.min,comparison geno.tau.R02<-get.tau(DengueSimR02,geno.type.func,r=r.max,r.low=r.min,comparison.type="independent") sero.tau.representative<-get.tau(DengueSimRepresentative, - sero.type.rep.func,r=r.max,r.low=r.min,comparison.type="representative") + sero.type.rep.func,r=r.max,r.low=r.min,comparison.type="representative") ## R0 of 1 plot(r.mid,sero.tau.R01,ylim=c(0.3,max(geno.tau.R01)),log="y", @@ -42,8 +42,8 @@ abline(v=100,lty=1,lwd=2) lines(r.mid,geno.tau.R01,pch=20,col=rgb(t(col2rgb("dark green")/255),alpha=0.6),lwd=1) lines(r.mid,sero.tau.representative,pch=20,col=rgb(t(col2rgb("dark blue")/255),alpha=0.6),lty=2) legend("topright",legend=c("Genotype","Serotype","Serotype (representative population)", - "Maximum transmission distance"),lwd=1,col=c("dark green","blue","blue","black"), - lty=c(1,1,2,1),bty="n") + "Maximum transmission distance"),lwd=1,col=c("dark green","blue","blue","black"), + lty=c(1,1,2,1),bty="n") ## R0 of 2 plot(r.mid,sero.tau.R02,ylim=c(0.3,max(geno.tau.R02)),log="y", diff --git a/R/spatialfuncs.r b/R/spatialfuncs.r index 6d0d4f7..21e5ea1 100644 --- a/R/spatialfuncs.r +++ b/R/spatialfuncs.r @@ -5,29 +5,32 @@ ##' specified by the passed in function with that point. ##' ##' @param posmat a matrix with columns x, y and any other named columns -##' columns needed by fun -##' @param fun a function that takes in two rows of posmat and returns: +##' columns needed by \code{fun} +##' @param fun a function that takes in two rows of \code{posmat} and returns: ##' \enumerate{ -##' \item for pairs included in the numerator and denominator +##' \item for pairs included in the numerator and denominator ##' \item for pairs that should only be included in the denominator ##' \item for pairs that should be ignored all together} ##' Note that names from \code{posmat} are not preserved in calls to \code{fun}, so the columns of the matrix should be ##' referenced numerically -##' so this is not available to the fun +##' so this is not available to the \code{fun} ##' @param r the series of spatial distances (or there maximums) we are ##' interested in ##' @param r.low the low end of each range, 0 by default ##' ##' @return pi value for each distance range that we look at. Where: ##' -##' \deqn{ \pi(d_1,d_2) = \frac{\sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j) \in \{1,2\}) }} +##'\deqn{ \pi(d_1, d_2) = \frac{\sum \boldsymbol{1} (d_{ij} \in [d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} [d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j) \in \{1,2\}) }} ##' ##' @author Justin Lessler and Henrik Salje ##' ##' @family get.pi ##' @family spatialtau ##' -##' @example R/examples/get_pi.R +##' @example +##' \dontrun{ +##' R/examples/get_pi.R +##' } ##' get.pi <- function(posmat, @@ -79,15 +82,19 @@ get.pi <- function(posmat, ##' ##' @return theta value for each distance range that we look at. Where: ##' -##' \deqn{ \theta(d_1,d_2) = \frac{\sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j)=2) }} +##' \deqn{ \theta(d_1,d_2) = \frac{\sum \boldsymbol{1} d_{ij} \in [d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} d_{ij} \in [d_1,d_2)) \boldsymbol{1} (f(i,j)=2) }} ##' ##' @author Justin Lessler and Henrik Salje ##' ##' @family get.theta ##' @family spatialtau ##' -##' @example R/examples/get_theta.R +##' @example +##' \dontrun{ +##' R/examples/get_theta.R +##' } ##' + get.theta <- function(posmat, fun, r = 1, @@ -132,7 +139,10 @@ get.theta <- function(posmat, ##' ##' @family get.pi ##' -##' @example R/examples/get_pi_typed.R +##' @example +##' \dontrun{ +##' R/examples/get_pi_typed.R +##' } ##' get.pi.typed <- function(posmat, typeA = -1, @@ -176,7 +186,10 @@ get.pi.typed <- function(posmat, ##' ##' @family get.theta ##' -##' @example R/examples/get_theta_typed.R +##' @example +##' \dontrun{ +##' R/examples/get_theta_typed.R +##' } ##' get.theta.typed <- function(posmat, typeA = -1, @@ -221,10 +234,10 @@ get.theta.typed <- function(posmat, ##' ##' @family get.pi ##' -##' @examples +##' @example ##' \dontrun{ -##' R/examples/get_pi_ci.R -##' } +##' R/examples/get_pi_ci.R +##' } ##' get.pi.ci <- function(posmat, fun, @@ -270,10 +283,10 @@ get.pi.ci <- function(posmat, ##' ##' @family get.theta ##' -##' @examples +##' @example ##' \dontrun{ -##' R/examples/get_theta_ci.R -##' } +##' R/examples/get_theta_ci.R +##' } ##' get.theta.ci <- function(posmat, fun, @@ -319,7 +332,7 @@ get.theta.ci <- function(posmat, ##' ##' @family get.pi ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_pi_bootstrap.R ##' } @@ -377,7 +390,7 @@ get.pi.bootstrap <- function(posmat, ##' ##' @family get.theta ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_theta_bootstrap.R ##' } @@ -430,7 +443,7 @@ get.theta.bootstrap <- function(posmat, ##' ##' @family get.pi ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_pi_typed_bootstrap.R ##' } @@ -481,7 +494,7 @@ get.pi.typed.bootstrap <- function(posmat, ##' ##' @family get.theta ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_theta_typed_bootstrap.R ##' } @@ -514,7 +527,6 @@ get.theta.typed.bootstrap <- function(posmat, return(rc) } - ##' get the null distribution of the \code{get.pi} function ##' ##' Does permutations to calculate the null distribution of get pi @@ -531,7 +543,7 @@ get.theta.typed.bootstrap <- function(posmat, ##' ##' @family get.pi ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_pi_permute.R ##' } @@ -586,7 +598,7 @@ get.pi.permute <- function(posmat, ##' ##' @family get.theta ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_theta_permute.R ##' } @@ -644,7 +656,7 @@ get.theta.permute <- function(posmat, ##' ##' @family get.pi ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_pi_typed_permute.R ##' } @@ -705,7 +717,7 @@ get.pi.typed.permute <- function(posmat, ##' ##' @family get.theta ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_theta_typed_permute.R ##' } @@ -776,18 +788,18 @@ get.theta.typed.permute <- function(posmat, ##' ##' @return The tau value for each distance we look at. If \code{comparison.type} is "representative", this is: ##' -##' \code{tau = get.pi(posmat, fun, r, r.low)/get.pi(posmat,fun,0,infinity)} +##' \code{tau = get.pi(posmat, fun, r, r.low)/get.pi(posmat,fun,infinity,0)} ##' ##' If \code{comparison.type} is "independent", this is: ##' -##' \code{tau = get.theta(posmat, fun, r, r.low)/get.theta(posmat,fun,0,infinity)} +##' \code{tau = get.theta(posmat, fun, r, r.low)/get.theta(posmat,fun,infinity,0)} ##' ##' @author Justin Lessler and Henrik Salje ##' ##' @family get.tau ##' @family spatialtau ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_tau.R ##' } @@ -827,7 +839,7 @@ get.tau <- function(posmat, return(rc) } ##' -##' Optimizewd version of \code{get.tau} for typed data +##' Optimized version of \code{get.tau} for typed data ##' ##' Version of th e \code{get.tau} function that is optimized for ##' statically typed data. That is data where we want the relationship between @@ -850,8 +862,10 @@ get.tau <- function(posmat, ##' ##' @family get.tau ##' -##' @example R/examples/get_tau_typed.R -##' +##' @example +##' \dontrun{ +##' R/examples/get_tau_typed.R +##'} get.tau.typed <- function(posmat, typeA = -1, typeB = -1, @@ -906,7 +920,7 @@ get.tau.typed <- function(posmat, ##' ##' @family get.tau ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_tau_ci.R ##' } @@ -955,7 +969,7 @@ get.tau.ci <- function(posmat, ##' ##' @family get.tau ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_tau_bootstrap.R ##' } @@ -1022,7 +1036,7 @@ get.tau.bootstrap <- function(posmat, ##' ##' @family get.tau ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_tau_typed_bootstrap.R ##' } @@ -1086,7 +1100,7 @@ get.tau.typed.bootstrap <- function(posmat, ##' ##' @family get.tau ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_tau_permute.R ##' } @@ -1159,7 +1173,7 @@ get.tau.permute <- function(posmat, ##' ##' @family get.tau ##' -##' @examples +##' @example ##' \dontrun{ ##' R/examples/get_tau_typed_permute.R ##' } diff --git a/R/transdistfuncs.r b/R/transdistfuncs.r index 2cbdf6d..8e64186 100644 --- a/R/transdistfuncs.r +++ b/R/transdistfuncs.r @@ -19,7 +19,10 @@ ##' ##' @author Justin Lessler, Henrik Salje, and John Giles ##' -##' @example R/examples/sim_epidemic.R +##' @example +##' \dontrun{ +##' R/examples/sim_epidemic.R +##' } ##' sim.epidemic <- function( @@ -120,9 +123,11 @@ sim.plot <- function(sim) { ##' ##' @family est.wt ##' -##' @example R/examples/est_wt_matrix.R +##' @example +##' \dontrun{ +##' R/examples/est_wt_matrix.R +##'} ##' - est.wt.matrix <- function( case.times, # a vector of case times (the time step in which each case occurs) gen.t.dist, # the generation time distribution @@ -152,7 +157,7 @@ est.wt.matrix <- function( ##' Estimate matrix of basic Wallinga-Teunis weights ##' ##' A function called by \code{est.wt.matrix}, which calculates the basic weights in the Wallinga-Teunis matrix given -##' the time of each case occurrence and the generation time distribution of the pathogen. Code adapted from the R0 package. +##' the time of each case occurrence and the generation time distribution of the pathogen. Code adapted from the \pkg{R0} package. ##' ##' @param case.times a vector giving the occurrence time for each case ##' @param gen.t.dist a vector giving the generation time distribution for the infecting pathogen @@ -168,7 +173,10 @@ est.wt.matrix <- function( ##' ##' @family est.wt ##' -##' @example R/examples/est_wt_matrix_weights.R +##' @example +##' \dontrun{ +##' R/examples/est_wt_matrix_weights.R +##'} ##' est.wt.matrix.weights <- function( @@ -224,8 +232,10 @@ est.wt.matrix.weights <- function( ##' ##' @family transdist ##' -##' @example R/examples/get_transdist_theta.R -##' +##' @example +##' \dontrun{ +##' R/examples/get_transdist_theta.R +##'} ##' # 'get.theta' is already taken @@ -313,7 +323,10 @@ get.transdist.theta <-function(wal.teun.mat, ##' ##' @family transdist ##' -##' @example R/examples/est_transdist_theta_weights.R +##' @example +##' \dontrun{ +##' R/examples/est_transdist_theta_weights.R +##'} ##' est.transdist.theta.weights <- function(case.times, @@ -379,7 +392,10 @@ est.transdist.theta.weights <- function(case.times, ##' @family est.wt ##' @family transdist ##' -##' @example R/examples/est_transdist.R +##' @example +##' \dontrun{ +##' R/examples/est_transdist.R +##'} ##' est.transdist <- function( @@ -523,7 +539,10 @@ est.transdist <- function( ##' ##' @family transdist ##' -##' @example R/examples/est_transdist_bootstrap_ci.R +##' @example +##' \dontrun{ +##' R/examples/est_transdist_bootstrap_ci.R +##' } ##' est.transdist.bootstrap.ci <- function( @@ -646,7 +665,10 @@ est.transdist.bootstrap.ci <- function( ##' ##' @family transdist ##' -##' @example R/examples/est_transdist_temporal.R +##' @example +##' \dontrun{ +##' R/examples/est_transdist_temporal.R +##'} ##' est.transdist.temporal <- function( @@ -770,7 +792,10 @@ est.transdist.temporal <- function( ##' ##' @family transdist ##' -##' @example R/examples/est_transdist_temporal_bootstrap_ci.R +##' @example +##' \dontrun{ +##' R/examples/est_transdist_temporal_bootstrap_ci.R +##'} ##' est.transdist.temporal.bootstrap.ci <- function( diff --git a/R/wrapperfuncs.r b/R/wrapperfuncs.r index afd585e..13f380a 100644 --- a/R/wrapperfuncs.r +++ b/R/wrapperfuncs.r @@ -1,6 +1,6 @@ ##' Cross type K function using homotypic and heterotypic case types ##' -##' A wrapper function of the \link[spatstat]{Kcross} function from the \pkg{spatstat} package that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type K-function based on user defined case type homology +##' A wrapper function of the \link[spatstat]{Kcross} function from the \pkg{spatstat} package (Baddeley et al. 2016) that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type K-function based on user defined case type homology ##' ##' @param epi.data a three-column numerical matrix that contains coordinates (\code{x} and \code{y}) for each case and information on case type (e.g. genotype or serotype). First two columns must be \code{x} and \code{y} ##' @param type an integer giving the column that contains information on case type. Must be an integer or a character @@ -9,7 +9,7 @@ ##' @param r a numeric vector giving the spatial distances ##' @param correction type of edge correction to be applied (default set to simple 'border' edge correction). See the \link[spatstat]{Kcross} function in the \pkg{spatstat} package for more details ##' -##' @return a data frame with a minimum of two columns giving the radius \code{r} and value of the K function, where the titel of the column the type of edge correction used +##' @return a data frame with a minimum of three columns giving the radius (\code{r}), the theoretical value of the K function for a Poisson process (\code{theo}), and value of the K function evaluated at radius \code{r}. The column name gives the type of edge correction used ##' ##' @author Justin Lessler, Henrik Salje, and John Giles ##' @@ -27,16 +27,16 @@ get.cross.K <- function(epi.data, # matrix containing xy coordinates and case ty correction='border' # type of edge correction to be applied ){ - if(is.matrix(epi.data) == FALSE && is.numeric(epi.data) == FALSE) stop('epi.data must be a numeric matrix') + if(is.matrix(epi.data) == FALSE | is.numeric(epi.data) == FALSE) stop('epi.data must be a numeric matrix') + tmp <- type > ncol(epi.data); if(tmp == TRUE) stop('type column defined is out of bounds') xy <- epi.data[, 1:2, drop=FALSE] case.types <- epi.data[, type] unique.case.types <- unique(case.types) if(is.null(het)) het <- unique.case.types[!(unique.case.types %in% hom)] - - hom <- as.numeric(hom) - het <- as.numeric(het) + if(all(hom %in% unique.case.types) == FALSE | all(het %in% unique.case.types) == FALSE) + stop('all homotypic and heterotypic case types are not found in case type data') message("Calculating cross type K function for:") message(paste(c("Homotypic case type(s):", hom, paste("(n=", length(which(case.types %in% hom)), ")", sep="")), @@ -48,24 +48,24 @@ get.cross.K <- function(epi.data, # matrix containing xy coordinates and case ty message(paste(c("Excluded case type(s):", excl, paste("(n=", length(which(case.types %in% excl)), ")", sep="")), collapse=" ")) - if(all(hom %in% unique.case.types) == FALSE | all(het %in% unique.case.types) == FALSE) - stop('all homotypic and heterotypic case types are not found in case type data') - case.marks <- rep(NA, length(case.types)) case.marks[which(case.types %in% het)] <- 0 case.marks[which(case.types %in% hom)] <- 1 ppp.dat <- cbind(xy, case.types, case.marks)[complete.cases(case.marks),] - epi.data.ppp <- spatstat::as.ppp(ppp.dat, spatstat::bounding.box.xy(xy), check=FALSE) - marks(epi.data.ppp) <- as.factor(ppp.dat[,ncol(ppp.dat)]) + epi.data.ppp <- spatstat::ppp(xy[,1], xy[,2], + spatstat::bounding.box.xy(xy), + marks=as.factor(ppp.dat[,ncol(ppp.dat)]), + check=FALSE) + #marks(epi.data.ppp) <- - return(as.data.frame(spatstat::Kcross(X=epi.data.ppp, i=0, j=1, r=r, correction=correction))[,-2]) + return(as.data.frame(spatstat::Kcross(X=epi.data.ppp, i=0, j=1, r=r, correction=correction))) } ##' Cross type Pair Correlation Function using homotypic and heterotypic case types ##' -##' A wrapper function of the \link[spatstat]{pcf} and \link[spatstat]{Kcross} functions from the \pkg{spatstat} package that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type Pair Correlation Function based on user defined case type homology +##' A wrapper function of the \link[spatstat]{pcf} and \link[spatstat]{Kcross} functions from the \pkg{spatstat} package (Baddeley et al. 2016) that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type Pair Correlation Function based on user defined case type homology ##' ##' @param epi.data a three-column numerical matrix that contains coordinates (\code{x} and \code{y}) for each case and information on case type (e.g. genotype or serotype). First two columns must be \code{x} and \code{y} ##' @param type an integer giving the column that contains information on case type. Must be an integer or a character @@ -74,7 +74,7 @@ get.cross.K <- function(epi.data, # matrix containing xy coordinates and case ty ##' @param r a numeric vector giving the spatial distances ##' @param correction type of edge correction to be applied (default set to simple 'border' edge correction). See the \link[spatstat]{Kcross} function in the \pkg{spatstat} package for more details ##' -##' @return a data frame with two columns giving the radius \code{r} and value of the Pair Correlation Function \code{pcf} +##' @return a data frame with two columns giving the radius \code{r}, the theoretical value of the Pair Correlation Function for a Poisson process (\code{theo}), and value of the Pair Correlation Function \code{pcf} ##' ##' @author Justin Lessler, Henrik Salje, and John Giles ##' @@ -92,18 +92,18 @@ get.cross.PCF <- function(epi.data, # matrix containing xy coordinates and case correction='border' # type of edge correction to be applied ){ - if(is.matrix(epi.data) == FALSE && is.numeric(epi.data) == FALSE) stop('epi.data must be a numeric matrix') + if(is.matrix(epi.data) == FALSE | is.numeric(epi.data) == FALSE) stop('epi.data must be a numeric matrix') + tmp <- type > ncol(epi.data); if(tmp == TRUE) stop('type column defined is out of bounds') xy <- epi.data[, 1:2, drop=FALSE] case.types <- epi.data[, type] unique.case.types <- unique(case.types) if(is.null(het)) het <- unique.case.types[!(unique.case.types %in% hom)] + if(all(hom %in% unique.case.types) == FALSE | all(het %in% unique.case.types) == FALSE) + stop('all homotypic and heterotypic case types are not found in case type data') - hom <- as.numeric(hom) - het <- as.numeric(het) - - message("Calculating cross type Pair Correlation Function for:") + message("Calculating cross type PCF function for:") message(paste(c("Homotypic case type(s):", hom, paste("(n=", length(which(case.types %in% hom)), ")", sep="")), collapse=" ")) message(paste(c("Heterotypic case type(s):", het, paste("(n=", length(which(case.types %in% het)), ")", sep="")), @@ -113,17 +113,17 @@ get.cross.PCF <- function(epi.data, # matrix containing xy coordinates and case message(paste(c("Excluded case type(s):", excl, paste("(n=", length(which(case.types %in% excl)), ")", sep="")), collapse=" ")) - if(all(hom %in% unique.case.types) == FALSE | all(het %in% unique.case.types) == FALSE) - stop('all homotypic and heterotypic case types are not found in case type data') - case.marks <- rep(NA, length(case.types)) case.marks[which(case.types %in% het)] <- 0 case.marks[which(case.types %in% hom)] <- 1 ppp.dat <- cbind(xy, case.types, case.marks)[complete.cases(case.marks),] - epi.data.ppp <- spatstat::as.ppp(ppp.dat, spatstat::bounding.box.xy(xy), check=FALSE) - marks(epi.data.ppp) <- as.factor(ppp.dat[,ncol(ppp.dat)]) + epi.data.ppp <- spatstat::ppp(xy[,1], xy[,2], + spatstat::bounding.box.xy(xy), + marks=as.factor(ppp.dat[,ncol(ppp.dat)]), + check=FALSE) - return(as.data.frame(spatstat::pcf(spatstat::Kcross(X=epi.data.ppp, i=0, j=1, r=r, correction=correction), spar=1, method='b'))[,-2]) -} \ No newline at end of file + return(as.data.frame(spatstat::pcf(spatstat::Kcross(X=epi.data.ppp, i=0, j=1, r=r, correction=correction), spar=0.8, method='c'))) +} + diff --git a/man/est.transdist.Rd b/man/est.transdist.Rd index 9fc2fb2..15e3710 100644 --- a/man/est.transdist.Rd +++ b/man/est.transdist.Rd @@ -32,30 +32,6 @@ a list containing the estimated mean distance of the transmission kernel (\code{ \description{ this function estimates the mean transmission distance of an epidemic when given the locations and times of symptomatic cases and the mean and standard deviation of the generation time of the infecting pathogen } -\examples{ -set.seed(123) - -# Exponentially distributed transmission kernel with mean and standard deviation = 100 -dist.func <- alist(n=1, a=1/100, rexp(n, a)) - -# Simulate epidemic -a <- sim.epidemic(R=1.5, - gen.t.mean=7, - gen.t.sd=2, - min.cases=50, - tot.generations=12, - trans.kern.func=dist.func) - -# Estimate mean and standara deviation of transmission kernel -b <- est.transdist(epi.data=a, - gen.t.mean=7, - gen.t.sd=2, - t1=0, - max.sep=1e10, - max.dist=1e10, - n.transtree.reps=10) -b -} \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. } @@ -72,3 +48,5 @@ Other transdist: \code{\link{est.transdist.bootstrap.ci}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{est.wt} +\concept{transdist} diff --git a/man/est.transdist.bootstrap.ci.Rd b/man/est.transdist.bootstrap.ci.Rd index 9272a29..834fcd5 100644 --- a/man/est.transdist.bootstrap.ci.Rd +++ b/man/est.transdist.bootstrap.ci.Rd @@ -43,35 +43,6 @@ a list object containing the point estimate for mean transmission distance and l } \description{ Runs \code{est.trandsdist} on multiple bootstraps of the data and calculates confidence intervals for the mean transmission distance. -} -\examples{ -set.seed(1) - -# Exponentially distributed transmission kernel with mean and standard deviation = 100 -dist.func <- alist(n=1, a=1/100, rexp(n, a)) - -# Simulate epidemic -a <- sim.epidemic(R=2.5, - gen.t.mean=7, - gen.t.sd=2, - min.cases=20, - tot.generations=5, - trans.kern.func=dist.func) - -# Estimate mean transmission kernel and its bootstrapped confidence intervals -b <- est.transdist.bootstrap.ci(epi.data=a, - gen.t.mean=7, - gen.t.sd=2, - t1=0, - max.sep=1e10, - max.dist=1e10, - n.transtree.reps=10, - mean.equals.sd=TRUE, - boot.iter=10, - ci.low=0.025, - ci.high=0.975, - parallel=FALSE) - } \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. @@ -86,3 +57,4 @@ Other transdist: \code{\link{est.transdist.temporal.bootstrap.ci}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{transdist} diff --git a/man/est.transdist.temporal.Rd b/man/est.transdist.temporal.Rd index bbba35d..004d93f 100644 --- a/man/est.transdist.temporal.Rd +++ b/man/est.transdist.temporal.Rd @@ -4,9 +4,9 @@ \alias{est.transdist.temporal} \title{Change in mean transmission distance over time} \usage{ -est.transdist.temporal(epi.data, gen.t.mean, gen.t.sd, t1, max.sep, max.dist, - n.transtree.reps = 10, mean.equals.sd = FALSE, theta.weights = NULL, - parallel = FALSE, n.cores = NULL) +est.transdist.temporal(epi.data, gen.t.mean, gen.t.sd, t1, max.sep, + max.dist, n.transtree.reps = 10, mean.equals.sd = FALSE, + theta.weights = NULL, parallel = FALSE, n.cores = NULL) } \arguments{ \item{epi.data}{a three-column matrix giving the coordinates (\code{x} and \code{y}) and time of infection (\code{t} for all cases in an epidemic (columns must be in \code{x}, \code{y}, \code{t} order)} @@ -39,42 +39,6 @@ NAs are returned for time steps which contain fewer than three cases Estimates the change in mean transmission distance over the duration of the epidemic by running \code{est.trandsdist} on all cases occuring up to each time point. } -\examples{ -set.seed(123) - -# Exponentially distributed transmission kernel with mean and standard deviation = 100 -dist.func <- alist(n=1, a=1/100, rexp(n, a)) - -# Simulate epidemic -a <- sim.epidemic(R=2, - gen.t.mean=7, - gen.t.sd=2, - tot.generations=7, - min.cases=30, - trans.kern.func=dist.func) - -a <- a[sample(1:nrow(a), 50),] # subsample a to 50 observations - -# Estimate mean transmission kernel over time -b <- est.transdist.temporal(epi.data=a, - gen.t.mean=7, - gen.t.sd=2, - t1=0, - max.sep=1e10, - max.dist=1e10, - n.transtree.reps=5, - mean.equals.sd=TRUE, - parallel=FALSE) - -plot(b[,1], pch=19, col='grey', ylim=c(min(b[,1], na.rm=TRUE), max(b[,1], na.rm=TRUE)), - xlab='Time step', ylab='Estimated mean of transmission kernel') -abline(h=100, col='red', lty=2) -axis(3, b[,2]) - -low <- loess(b[,1] ~ as.vector(1:length(b[,1]))) -low <- predict(low, newdata=data.frame(as.vector(1:length(b[,1])))) -lines(low, lwd=3, col='blue') -} \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. } @@ -88,3 +52,4 @@ Other transdist: \code{\link{est.transdist.bootstrap.ci}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{transdist} diff --git a/man/est.transdist.temporal.bootstrap.ci.Rd b/man/est.transdist.temporal.bootstrap.ci.Rd index 132c862..9ecdbb6 100644 --- a/man/est.transdist.temporal.bootstrap.ci.Rd +++ b/man/est.transdist.temporal.bootstrap.ci.Rd @@ -4,8 +4,8 @@ \alias{est.transdist.temporal.bootstrap.ci} \title{Bootstrapped confidence intervals for the change in mean transmission distance over time} \usage{ -est.transdist.temporal.bootstrap.ci(epi.data, gen.t.mean, gen.t.sd, t1, max.sep, - max.dist, n.transtree.reps = 100, mean.equals.sd = FALSE, +est.transdist.temporal.bootstrap.ci(epi.data, gen.t.mean, gen.t.sd, t1, + max.sep, max.dist, n.transtree.reps = 100, mean.equals.sd = FALSE, theta.weights = NULL, boot.iter, ci.low = 0.025, ci.high = 0.975, parallel = FALSE, n.cores = NULL) } @@ -45,56 +45,6 @@ a four-column numeric matrix containing the point estimate for mean transmission Estimates bootstrapped confidence intervals for the mean transmission distance over the duration of the epidemic by running \code{est.trandsdist} on all cases occuring up to each time point. } -\examples{ -\donttest{ - set.seed(123) - - # Exponentially distributed transmission kernel with mean and standard deviation = 100 - dist.func <- alist(n=1, a=1/100, rexp(n, a)) - - # Simulate epidemic - a <- sim.epidemic(R=2, - gen.t.mean=7, - gen.t.sd=2, - tot.generations=8, - min.cases=30, - trans.kern.func=dist.func) - - a <- a[sample(1:nrow(a), 100),] # subsample a to 100 observations - - # Estimate change in mean transmission kernel over time with confidence intervals - nc <- 4 # Run in parallel - - b <- est.transdist.temporal.bootstrap.ci(epi.data=a, - gen.t.mean=7, - gen.t.sd=2, - t1=0, - max.sep=1e10, - max.dist=1e10, - n.transtree.reps=10, - mean.equals.sd=TRUE, - boot.iter=10, - ci.low=0.025, - ci.high=0.975, - parallel=TRUE, - n.cores=nc) - - plot(b[,1], pch=19, col='grey', ylim=c(min(b[,1:3], na.rm=TRUE), max(b[,1:3], na.rm=TRUE)), - xlab='Time step', ylab='Estimated mean of transmission kernel') - abline(h=100, col='red', lty=2) - axis(3, 1:nrow(b), b[,4]) - - low <- loess(b[,1] ~ as.vector(1:nrow(b)), span=1) - low <- predict(low, newdata=data.frame(as.vector(1:nrow(b)))) - lines(low, lwd=3, col='blue') - - for(i in 2:3) { - low <- loess(b[,i] ~ as.vector(1:nrow(b)), span=1) - low <- predict(low, newdata=data.frame(as.vector(1:nrow(b)))) - lines(low, lty=2, lwd=3, col='blue') - } -} -} \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. } @@ -108,3 +58,4 @@ Other transdist: \code{\link{est.transdist.bootstrap.ci}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{transdist} diff --git a/man/est.transdist.theta.weights.Rd b/man/est.transdist.theta.weights.Rd index cd6cd79..07e0377 100644 --- a/man/est.transdist.theta.weights.Rd +++ b/man/est.transdist.theta.weights.Rd @@ -25,46 +25,6 @@ a three-dimensional array containing the mean normalized theta weights estimated This function estimates the weight of each theta value by performing a user defined number of replications with the \code{get.transdist.theta} function. The weights of each theta are calculated as the number of simulations in which a case at time \code{t1} and \code{t2} are separated by theta transmission events. } -\examples{ -set.seed(1) - -gen.t.mean <- 7 -gen.t.sd <- 2 -t1 <- 0 - -# Normally distributed transmission kernel with mean and standard deviation = 100 -dist.func <- alist(n=1, a=1/100, rexp(n, a)) - -# Simulate epidemic -a <- sim.epidemic(R=5, - gen.t.mean=gen.t.mean, - gen.t.sd=gen.t.sd, - min.cases=5, - tot.generations=3, - trans.kern.func=dist.func) - -# Get case times -a <- a[order(a[,3]),] -case.times <- round(a[,3]) -unique.times <- unique(case.times) -ntimes <- length(unique.times) - - -# Generation time distribution -max.t <- round(max(unique.times) - t1) - 1 -n.step <- round(max.t/gen.t.mean) -gen <- rep(0, max.t*2) -for (i in 1:n.step){gen <- gen + dnorm(1:(max.t*2), gen.t.mean*i, gen.t.sd*i)} -gen[1] <- 0 # No instantaneous infections -t.density <- gen/sum(gen) - -# Estimation of theta weights matrix -b <- est.transdist.theta.weights(case.times=case.times, - n.rep=3, - gen.t.mean=gen.t.mean, - t1=t1, - t.density=t.density) -} \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. } @@ -78,3 +38,4 @@ Other transdist: \code{\link{est.transdist.bootstrap.ci}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{transdist} diff --git a/man/est.wt.matrix.Rd b/man/est.wt.matrix.Rd index 30d91d1..95fd7c5 100644 --- a/man/est.wt.matrix.Rd +++ b/man/est.wt.matrix.Rd @@ -22,13 +22,6 @@ A function which takes the time of each case occurrence, the generation time dis and the matrix of basic Wallinga-Teunis weights and estimates the probability that an infectee occurring time step j (columns) was infected by a case occurring at time i (rows). } -\examples{ -case.times <- c(1,2,2,3,3) -gen <- c(0, 2/3, 1/3, 0, 0) -t.density <- gen/sum(gen) - -a <- est.wt.matrix(case.times=case.times, gen.t.dist=t.density) -} \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. } @@ -39,3 +32,4 @@ Other est.wt: \code{\link{est.transdist}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{est.wt} diff --git a/man/est.wt.matrix.weights.Rd b/man/est.wt.matrix.weights.Rd index 5975804..205ff2a 100644 --- a/man/est.wt.matrix.weights.Rd +++ b/man/est.wt.matrix.weights.Rd @@ -16,14 +16,7 @@ a numerical matrix with the number of columns and rows equal to the number of ti } \description{ A function called by \code{est.wt.matrix}, which calculates the basic weights in the Wallinga-Teunis matrix given -the time of each case occurrence and the generation time distribution of the pathogen. Code adapted from the R0 package. -} -\examples{ -case.times <- c(1,2,2,3,3) -gen <- c(0, 2/3, 1/3, 0, 0) -t.density <- gen/sum(gen) - -a <- est.wt.matrix.weights(case.times=case.times, gen.t.dist=t.density) +the time of each case occurrence and the generation time distribution of the pathogen. Code adapted from the \pkg{R0} package. } \references{ Boelle P and Obadia T (2015). R0: Estimation of R0 and Real-Time Reproduction Number from Epidemics. R package version 1.2-6, \href{https://CRAN.R-project.org/package=R0}{https://CRAN.R-project.org/package=R0}. @@ -37,3 +30,4 @@ Other est.wt: \code{\link{est.transdist}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{est.wt} diff --git a/man/get.cross.K.Rd b/man/get.cross.K.Rd index db1b967..8dc1453 100644 --- a/man/get.cross.K.Rd +++ b/man/get.cross.K.Rd @@ -21,18 +21,17 @@ get.cross.K(epi.data, type, hom, het = NULL, r = NULL, \item{correction}{type of edge correction to be applied (default set to simple 'border' edge correction). See the \link[spatstat]{Kcross} function in the \pkg{spatstat} package for more details} } \value{ -a data frame with a minimum of two columns giving the radius \code{r} and value of the K function, where the titel of the column the type of edge correction used +a data frame with a minimum of three columns giving the radius (\code{r}), the theoretical value of the K function for a Poisson process (\code{theo}), and value of the K function evaluated at radius \code{r}. The column name gives the type of edge correction used } \description{ -A wrapper function of the \link[spatstat]{Kcross} function from the \pkg{spatstat} package that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type K-function based on user defined case type homology +A wrapper function of the \link[spatstat]{Kcross} function from the \pkg{spatstat} package (Baddeley et al. 2016) that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type K-function based on user defined case type homology } \examples{ data(DengueSimulationR01) k <- get.cross.K(epi.data=DengueSimR01, type=5, hom=2, het=NULL, r=NULL, correction='border') -k.theo <- pi*(k$r)^2 -plot(k.theo, type='l', col='red', lty=2, xlab='r', ylab='cross K function') +plot(k[,2], type='l', col='red', lty=2, xlab='r', ylab='cross K function') lines(k$border) } \references{ diff --git a/man/get.cross.PCF.Rd b/man/get.cross.PCF.Rd index fb9e336..25abb98 100644 --- a/man/get.cross.PCF.Rd +++ b/man/get.cross.PCF.Rd @@ -21,10 +21,10 @@ get.cross.PCF(epi.data, type, hom, het = NULL, r = NULL, \item{correction}{type of edge correction to be applied (default set to simple 'border' edge correction). See the \link[spatstat]{Kcross} function in the \pkg{spatstat} package for more details} } \value{ -a data frame with two columns giving the radius \code{r} and value of the Pair Correlation Function \code{pcf} +a data frame with two columns giving the radius \code{r}, the theoretical value of the Pair Correlation Function for a Poisson process (\code{theo}), and value of the Pair Correlation Function \code{pcf} } \description{ -A wrapper function of the \link[spatstat]{pcf} and \link[spatstat]{Kcross} functions from the \pkg{spatstat} package that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type Pair Correlation Function based on user defined case type homology +A wrapper function of the \link[spatstat]{pcf} and \link[spatstat]{Kcross} functions from the \pkg{spatstat} package (Baddeley et al. 2016) that takes epidemiological data used by \pkg{IDSpatialStats} functions and calculates the cross type Pair Correlation Function based on user defined case type homology } \examples{ data(DengueSimulationR01) diff --git a/man/get.pi.Rd b/man/get.pi.Rd index 7591ec5..8897ad5 100644 --- a/man/get.pi.Rd +++ b/man/get.pi.Rd @@ -8,16 +8,16 @@ get.pi(posmat, fun, r = 1, r.low = rep(0, length(r))) } \arguments{ \item{posmat}{a matrix with columns x, y and any other named columns -columns needed by fun} +columns needed by \code{fun}} -\item{fun}{a function that takes in two rows of posmat and returns: +\item{fun}{a function that takes in two rows of \code{posmat} and returns: \enumerate{ - \item for pairs included in the numerator and denominator + \item for pairs included in the numerator and denominator \item for pairs that should only be included in the denominator \item for pairs that should be ignored all together} Note that names from \code{posmat} are not preserved in calls to \code{fun}, so the columns of the matrix should be referenced numerically -so this is not available to the fun} +so this is not available to the \code{fun}} \item{r}{the series of spatial distances (or there maximums) we are interested in} @@ -27,27 +27,13 @@ interested in} \value{ pi value for each distance range that we look at. Where: -\deqn{ \pi(d_1,d_2) = \frac{\sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j) \in \{1,2\}) }} +\deqn{ \pi(d_1, d_2) = \frac{\sum \boldsymbol{1} (d_{ij} \in [d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} [d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j) \in \{1,2\}) }} } \description{ Generalized version of the \code{get.pi} function that takes in an arbitrary function and returns the probability that a point within a particular range of a point of interest shares the relationship specified by the passed in function with that point. } -\examples{ -data(DengueSimulationR02) - -r.max<-seq(20,1000,20) -r.min<-seq(0,980,20) - -sero.type.func<-function(a,b,tlimit=20){ - if(a[5]==b[5]&(abs(a[3]-b[3])<=tlimit)){rc=1} - else{rc=2} - return(rc) -} - -sero.pi<-get.pi(DengueSimR02,sero.type.func,r=r.max,r.low=r.min) -} \seealso{ Other get.pi: \code{\link{get.pi.bootstrap}}, \code{\link{get.pi.ci}}, \code{\link{get.pi.permute}}, @@ -61,3 +47,5 @@ Other spatialtau: \code{\link{get.tau}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.pi} +\concept{spatialtau} diff --git a/man/get.pi.bootstrap.Rd b/man/get.pi.bootstrap.Rd index af91509..1ec8a70 100644 --- a/man/get.pi.bootstrap.Rd +++ b/man/get.pi.bootstrap.Rd @@ -30,12 +30,6 @@ points and themselves will not be calculated. In each bootstrap iteration N observations are drawn from the existing data with replacement. To avoid errors in inference resulting from the same observatin being compared with itself in the bootstrapped data set, original indices are perserved, and pairs of points in the bootstrapped dataset with the same original index are ignored. -} -\examples{ -\dontrun{ - R/examples/get_pi_bootstrap.R - } - } \seealso{ Other get.pi: \code{\link{get.pi.ci}}, @@ -47,3 +41,4 @@ Other get.pi: \code{\link{get.pi.ci}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.pi} diff --git a/man/get.pi.ci.Rd b/man/get.pi.ci.Rd index 66ceb8d..23c1be4 100644 --- a/man/get.pi.ci.Rd +++ b/man/get.pi.ci.Rd @@ -4,8 +4,8 @@ \alias{get.pi.ci} \title{Calculate bootstrapped confidence intervals for \code{get.pi} values.} \usage{ -get.pi.ci(posmat, fun, r = 1, r.low = rep(0, length(r)), boot.iter = 1000, - ci.low = 0.025, ci.high = 0.975) +get.pi.ci(posmat, fun, r = 1, r.low = rep(0, length(r)), + boot.iter = 1000, ci.low = 0.025, ci.high = 0.975) } \arguments{ \item{posmat}{a matrix with columns type, x and y} @@ -29,12 +29,6 @@ a matrix with a row for the high and low values and \description{ Wrapper to \code{get.pi.bootstrap} that takes care of calculating the confidence intervals based on the bootstrapped values.. -} -\examples{ -\dontrun{ - R/examples/get_pi_ci.R - } - } \seealso{ Other get.pi: \code{\link{get.pi.bootstrap}}, @@ -46,3 +40,4 @@ Other get.pi: \code{\link{get.pi.bootstrap}}, \author{ Justin Lessler } +\concept{get.pi} diff --git a/man/get.pi.permute.Rd b/man/get.pi.permute.Rd index 8e6c32a..e13b63b 100644 --- a/man/get.pi.permute.Rd +++ b/man/get.pi.permute.Rd @@ -4,7 +4,8 @@ \alias{get.pi.permute} \title{get the null distribution of the \code{get.pi} function} \usage{ -get.pi.permute(posmat, fun, r = 1, r.low = rep(0, length(r)), permutations) +get.pi.permute(posmat, fun, r = 1, r.low = rep(0, length(r)), + permutations) } \arguments{ \item{posmat}{a matrix with columns type, x and y} @@ -24,12 +25,6 @@ pi values for all the distances we looked at Does permutations to calculate the null distribution of get pi if there were no spatial dependence. Randomly reassigns coordinates to each observation permutations times -} -\examples{ -\dontrun{ - R/examples/get_pi_permute.R - } - } \seealso{ Other get.pi: \code{\link{get.pi.bootstrap}}, @@ -38,3 +33,4 @@ Other get.pi: \code{\link{get.pi.bootstrap}}, \code{\link{get.pi.typed.permute}}, \code{\link{get.pi.typed}}, \code{\link{get.pi}} } +\concept{get.pi} diff --git a/man/get.pi.typed.Rd b/man/get.pi.typed.Rd index fbea488..01b913d 100644 --- a/man/get.pi.typed.Rd +++ b/man/get.pi.typed.Rd @@ -26,18 +26,6 @@ Version of the \code{get.pi} function that is optimized for statically typed dat data where we are interested in the probability of points within some distance of points of typeA are of typeB. } -\examples{ -data(DengueSimulationR02) - -r.max<-seq(20,1000,20) -r.min<-seq(0,980,20) - -#Lets see if there's a difference in spatial dependence by time case occurs -type<-2-(DengueSimR02[,"time"]<120) -tmp<-cbind(DengueSimR02,type=type) - -typed.pi<-get.pi.typed(tmp,typeA=1,typeB=2,r=r.max,r.low=r.min) -} \seealso{ Other get.pi: \code{\link{get.pi.bootstrap}}, \code{\link{get.pi.ci}}, \code{\link{get.pi.permute}}, @@ -47,3 +35,4 @@ Other get.pi: \code{\link{get.pi.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.pi} diff --git a/man/get.pi.typed.bootstrap.Rd b/man/get.pi.typed.bootstrap.Rd index d3fd70f..83ebd5b 100644 --- a/man/get.pi.typed.bootstrap.Rd +++ b/man/get.pi.typed.bootstrap.Rd @@ -26,12 +26,6 @@ pi values for all the distances we looked at \description{ Bootstraps typed pi values. Makes sure distances between a sample and another draw of itself are left out -} -\examples{ -\dontrun{ - R/examples/get_pi_typed_bootstrap.R - } - } \seealso{ Other get.pi: \code{\link{get.pi.bootstrap}}, @@ -39,3 +33,4 @@ Other get.pi: \code{\link{get.pi.bootstrap}}, \code{\link{get.pi.typed.permute}}, \code{\link{get.pi.typed}}, \code{\link{get.pi}} } +\concept{get.pi} diff --git a/man/get.pi.typed.permute.Rd b/man/get.pi.typed.permute.Rd index 3705b19..05a3dda 100644 --- a/man/get.pi.typed.permute.Rd +++ b/man/get.pi.typed.permute.Rd @@ -27,12 +27,6 @@ pi values for all the distances we looked at Does permutations to calculate the null distribution of get pi if there were no spatial dependence. Randomly reassigns coordinates to each observation permutations times -} -\examples{ -\dontrun{ - R/examples/get_pi_typed_permute.R - } - } \seealso{ Other get.pi: \code{\link{get.pi.bootstrap}}, @@ -43,3 +37,4 @@ Other get.pi: \code{\link{get.pi.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.pi} diff --git a/man/get.tau.Rd b/man/get.tau.Rd index d596408..b7bc5a7 100644 --- a/man/get.tau.Rd +++ b/man/get.tau.Rd @@ -34,22 +34,16 @@ interested in} \value{ The tau value for each distance we look at. If \code{comparison.type} is "representative", this is: -\code{tau = get.pi(posmat, fun, r, r.low)/get.pi(posmat,fun,0,infinity)} +\code{tau = get.pi(posmat, fun, r, r.low)/get.pi(posmat,fun,infinity,0)} If \code{comparison.type} is "independent", this is: -\code{tau = get.theta(posmat, fun, r, r.low)/get.theta(posmat,fun,0,infinity)} +\code{tau = get.theta(posmat, fun, r, r.low)/get.theta(posmat,fun,infinity,0)} } \description{ returns the relative probability (or odds) that points at some distance from an index point share some relationship with that point versus the probability (or odds) any point shares that relationship with that point. -} -\examples{ -\dontrun{ -R/examples/get_tau.R -} - } \seealso{ Other get.tau: \code{\link{get.tau.bootstrap}}, @@ -64,3 +58,5 @@ Other spatialtau: \code{\link{get.pi}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} +\concept{spatialtau} diff --git a/man/get.tau.bootstrap.Rd b/man/get.tau.bootstrap.Rd index 0045cc9..312febe 100644 --- a/man/get.tau.bootstrap.Rd +++ b/man/get.tau.bootstrap.Rd @@ -4,8 +4,8 @@ \alias{get.tau.bootstrap} \title{Bootstrap \code{get.tau} values.} \usage{ -get.tau.bootstrap(posmat, fun, r = 1, r.low = rep(0, length(r)), boot.iter, - comparison.type = "representative") +get.tau.bootstrap(posmat, fun, r = 1, r.low = rep(0, length(r)), + boot.iter, comparison.type = "representative") } \arguments{ \item{posmat}{a matrix appropriate for input to \code{get.tau}} @@ -27,12 +27,6 @@ tau values for all the distances we looked at Runs \code{get.tau} on multiple bootstraps of the data. Is formulated such that the relationship between points and themselves will not be calculated -} -\examples{ -\dontrun{ - R/examples/get_tau_bootstrap.R - } - } \seealso{ Other get.tau: \code{\link{get.tau.ci}}, @@ -44,3 +38,4 @@ Other get.tau: \code{\link{get.tau.ci}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} diff --git a/man/get.tau.ci.Rd b/man/get.tau.ci.Rd index 30665fb..987cb32 100644 --- a/man/get.tau.ci.Rd +++ b/man/get.tau.ci.Rd @@ -5,8 +5,8 @@ \title{Bootstrap confidence interval for the \code{get.tau} values} \usage{ get.tau.ci(posmat, fun, r = 1, r.low = rep(0, length(r)), - boot.iter = 1000, comparison.type = "representative", ci.low = 0.025, - ci.high = 0.975) + boot.iter = 1000, comparison.type = "representative", + ci.low = 0.025, ci.high = 0.975) } \arguments{ \item{posmat}{a matrix appropriate for input to \code{get.tau}} @@ -31,12 +31,6 @@ tau values for all the distances examined \description{ Wrapper to \code{get.tau.bootstrap} that takes care of calulating the confidence intervals based on the bootstrapped values -} -\examples{ -\dontrun{ - R/examples/get_tau_ci.R - } - } \seealso{ Other get.tau: \code{\link{get.tau.bootstrap}}, @@ -48,3 +42,4 @@ Other get.tau: \code{\link{get.tau.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} diff --git a/man/get.tau.permute.Rd b/man/get.tau.permute.Rd index 4537461..48c208f 100644 --- a/man/get.tau.permute.Rd +++ b/man/get.tau.permute.Rd @@ -4,8 +4,8 @@ \alias{get.tau.permute} \title{get the null distribution of the \code{get.tau} function} \usage{ -get.tau.permute(posmat, fun, r = 1, r.low = rep(0, length(r)), permutations, - comparison.type = "representative") +get.tau.permute(posmat, fun, r = 1, r.low = rep(0, length(r)), + permutations, comparison.type = "representative") } \arguments{ \item{posmat}{a matrix appropriate for input to \code{get.tau}} @@ -27,12 +27,6 @@ tau values for all the distances we looked at Does permutations to calculate the null distribution of get pi if there were no spatial dependence. Randomly reassigns coordinates to each observation permutations times -} -\examples{ -\dontrun{ - R/examples/get_tau_permute.R - } - } \seealso{ Other get.tau: \code{\link{get.tau.bootstrap}}, @@ -44,3 +38,4 @@ Other get.tau: \code{\link{get.tau.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} diff --git a/man/get.tau.typed.Rd b/man/get.tau.typed.Rd index 8373512..b637870 100644 --- a/man/get.tau.typed.Rd +++ b/man/get.tau.typed.Rd @@ -2,7 +2,7 @@ % Please edit documentation in R/spatialfuncs.r \name{get.tau.typed} \alias{get.tau.typed} -\title{Optimizewd version of \code{get.tau} for typed data} +\title{Optimized version of \code{get.tau} for typed data} \usage{ get.tau.typed(posmat, typeA = -1, typeB = -1, r = 1, r.low = rep(0, length(r)), comparison.type = "representative") @@ -32,23 +32,6 @@ Version of th e \code{get.tau} function that is optimized for statically typed data. That is data where we want the relationship between points of type A and points of type B } -\examples{ -data(DengueSimulationR02) - -r.max<-seq(20,1000,20) -r.min<-seq(0,980,20) -r.mid<-(r.max+r.min)/2 - -#Lets see if there's a difference in spatial dependence by time case occurs -type<-2-(DengueSimR02[,"time"]<120) -tmp<-cbind(DengueSimR02,type=type) - -typed.tau<-get.tau.typed(tmp,typeA=1,typeB=2,r=r.max,r.low=r.min,comparison.type = "independent") - -plot(r.mid,typed.tau,log="y",cex.axis=1.25, - xlab="Distance (m)",ylab="Tau",cex.main=0.9,lwd=2,type="l") -abline(h=1,lty=2) -} \seealso{ Other get.tau: \code{\link{get.tau.bootstrap}}, \code{\link{get.tau.ci}}, \code{\link{get.tau.permute}}, @@ -59,3 +42,4 @@ Other get.tau: \code{\link{get.tau.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} diff --git a/man/get.tau.typed.bootstrap.Rd b/man/get.tau.typed.bootstrap.Rd index 36c7162..8fd3957 100644 --- a/man/get.tau.typed.bootstrap.Rd +++ b/man/get.tau.typed.bootstrap.Rd @@ -5,7 +5,8 @@ \title{runs bootstrapping for \code{get.tau.typed}} \usage{ get.tau.typed.bootstrap(posmat, typeA = -1, typeB = -1, r = 1, - r.low = rep(0, length(r)), boot.iter, comparison.type = "representative") + r.low = rep(0, length(r)), boot.iter, + comparison.type = "representative") } \arguments{ \item{posmat}{a matrix with columns type, x and y} @@ -31,12 +32,6 @@ tau values for all the distances we looked at } \description{ runs bootstrapping for \code{get.tau.typed} -} -\examples{ -\dontrun{ - R/examples/get_tau_typed_bootstrap.R - } - } \seealso{ Other get.tau: \code{\link{get.tau.bootstrap}}, @@ -47,3 +42,4 @@ Other get.tau: \code{\link{get.tau.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} diff --git a/man/get.tau.typed.permute.Rd b/man/get.tau.typed.permute.Rd index ebe3968..493bdc5 100644 --- a/man/get.tau.typed.permute.Rd +++ b/man/get.tau.typed.permute.Rd @@ -32,12 +32,6 @@ a matrix with permutation tau values for each distance specified } \description{ get the null distribution for the \code{get.tau.typed} function -} -\examples{ -\dontrun{ -R/examples/get_tau_typed_permute.R -} - } \seealso{ Other get.tau: \code{\link{get.tau.bootstrap}}, @@ -48,3 +42,4 @@ Other get.tau: \code{\link{get.tau.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.tau} diff --git a/man/get.theta.Rd b/man/get.theta.Rd index 7f6bf11..f9cdf13 100644 --- a/man/get.theta.Rd +++ b/man/get.theta.Rd @@ -27,27 +27,13 @@ interested in} \value{ theta value for each distance range that we look at. Where: -\deqn{ \theta(d_1,d_2) = \frac{\sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} (d_{ij} \in (d_1,d_2)) \boldsymbol{1} (f(i,j)=2) }} +\deqn{ \theta(d_1,d_2) = \frac{\sum \boldsymbol{1} d_{ij} \in [d_1,d_2)) \boldsymbol{1} (f(i,j)=1) }{\sum \sum \boldsymbol{1} d_{ij} \in [d_1,d_2)) \boldsymbol{1} (f(i,j)=2) }} } \description{ Generalized version of the \code{get.theta} function that takes in an arbitrary function and returns the odds that a point within a particular range of a point of interest shares the relationship specified by the passed in function with that point. } -\examples{ -data(DengueSimulationR02) - -r.max<-seq(20,1000,20) -r.min<-seq(0,980,20) - -sero.type.func<-function(a,b,tlimit=20){ - if(a[5]==b[5]&(abs(a[3]-b[3])<=tlimit)){rc=1} - else{rc=2} - return(rc) -} - -sero.theta<-get.theta(DengueSimR02,sero.type.func,r=r.max,r.low=r.min) -} \seealso{ Other get.theta: \code{\link{get.theta.bootstrap}}, \code{\link{get.theta.ci}}, @@ -62,3 +48,5 @@ Other spatialtau: \code{\link{get.pi}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.theta} +\concept{spatialtau} diff --git a/man/get.theta.bootstrap.Rd b/man/get.theta.bootstrap.Rd index ad6455b..a453db4 100644 --- a/man/get.theta.bootstrap.Rd +++ b/man/get.theta.bootstrap.Rd @@ -30,12 +30,6 @@ points and themselves will not be calculated. In each bootstrap iteration N observations are drawn from the existing data with replacement. To avoid errors in inference resulting from the same observatin being compared with itself in the bootstrapped data set, original indices are perserved, and pairs of points in the bootstrapped dataset with the same original index are ignored. -} -\examples{ -\dontrun{ - R/examples/get_theta_bootstrap.R - } - } \seealso{ Other get.theta: \code{\link{get.theta.ci}}, @@ -47,3 +41,4 @@ Other get.theta: \code{\link{get.theta.ci}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.theta} diff --git a/man/get.theta.ci.Rd b/man/get.theta.ci.Rd index c9706d6..ebdd43b 100644 --- a/man/get.theta.ci.Rd +++ b/man/get.theta.ci.Rd @@ -29,12 +29,6 @@ a matrix with a row for the high and low values and \description{ Wrapper to \code{get.theta.bootstrap} that takes care of calculating the confience intervals based on the bootstrapped values. -} -\examples{ -\dontrun{ - R/examples/get_theta_ci.R - } - } \seealso{ Other get.theta: \code{\link{get.theta.bootstrap}}, @@ -46,3 +40,4 @@ Other get.theta: \code{\link{get.theta.bootstrap}}, \author{ Justin Lessler } +\concept{get.theta} diff --git a/man/get.theta.permute.Rd b/man/get.theta.permute.Rd index 80577b4..a67fb01 100644 --- a/man/get.theta.permute.Rd +++ b/man/get.theta.permute.Rd @@ -25,12 +25,6 @@ theta values for all the distances we looked at Does permutations to calculate the null distribution of get theta if there were no spatial dependence. Randomly reassigns coordinates to each observation permutations times -} -\examples{ -\dontrun{ - R/examples/get_theta_permute.R - } - } \seealso{ Other get.theta: \code{\link{get.theta.bootstrap}}, @@ -39,3 +33,4 @@ Other get.theta: \code{\link{get.theta.bootstrap}}, \code{\link{get.theta.typed.permute}}, \code{\link{get.theta.typed}}, \code{\link{get.theta}} } +\concept{get.theta} diff --git a/man/get.theta.typed.Rd b/man/get.theta.typed.Rd index 3cfce27..95d7b93 100644 --- a/man/get.theta.typed.Rd +++ b/man/get.theta.typed.Rd @@ -4,8 +4,8 @@ \alias{get.theta.typed} \title{Optimized version of \code{get.theta} for typed data.} \usage{ -get.theta.typed(posmat, typeA = -1, typeB = -1, r = 1, r.low = rep(0, - length(r))) +get.theta.typed(posmat, typeA = -1, typeB = -1, r = 1, + r.low = rep(0, length(r))) } \arguments{ \item{posmat}{a matrix with columns type, x and y} @@ -26,18 +26,6 @@ Version of the \code{get.theta} function that is optimized for statically typed data where we are interested in the odds that points within some distance of points of typeA are of typeB. } -\examples{ -data(DengueSimulationR02) - -r.max<-seq(20,1000,20) -r.min<-seq(0,980,20) - -#Lets see if there's a difference in spatial dependence by time case occurs -type<-2-(DengueSimR02[,"time"]<120) -tmp<-cbind(DengueSimR02,type=type) - -typed.theta.R01<-get.theta.typed(tmp,typeA=2,typeB=2,r=r.max,r.low=r.min) -} \seealso{ Other get.theta: \code{\link{get.theta.bootstrap}}, \code{\link{get.theta.ci}}, @@ -49,3 +37,4 @@ Other get.theta: \code{\link{get.theta.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.theta} diff --git a/man/get.theta.typed.bootstrap.Rd b/man/get.theta.typed.bootstrap.Rd index 7634087..eed3174 100644 --- a/man/get.theta.typed.bootstrap.Rd +++ b/man/get.theta.typed.bootstrap.Rd @@ -26,12 +26,6 @@ theta values for all the distances we looked at \description{ Bootstraps typed pi values. Makes sure distances between a sample and another draw of itself are left out -} -\examples{ -\dontrun{ - R/examples/get_theta_typed_bootstrap.R - } - } \seealso{ Other get.theta: \code{\link{get.theta.bootstrap}}, @@ -40,3 +34,4 @@ Other get.theta: \code{\link{get.theta.bootstrap}}, \code{\link{get.theta.typed.permute}}, \code{\link{get.theta.typed}}, \code{\link{get.theta}} } +\concept{get.theta} diff --git a/man/get.theta.typed.permute.Rd b/man/get.theta.typed.permute.Rd index 9d36168..e645bc7 100644 --- a/man/get.theta.typed.permute.Rd +++ b/man/get.theta.typed.permute.Rd @@ -27,12 +27,6 @@ theta values for all the distances we looked at Does permutations to calculate the null distribution of get theta if there were no spatial dependence. Randomly reassigns coordinates to each observation permutations times -} -\examples{ -\dontrun{ - R/examples/get_theta_typed_permute.R - } - } \seealso{ Other get.theta: \code{\link{get.theta.bootstrap}}, @@ -44,3 +38,4 @@ Other get.theta: \code{\link{get.theta.bootstrap}}, \author{ Justin Lessler and Henrik Salje } +\concept{get.theta} diff --git a/man/get.transdist.theta.Rd b/man/get.transdist.theta.Rd index 576fab6..926a119 100644 --- a/man/get.transdist.theta.Rd +++ b/man/get.transdist.theta.Rd @@ -25,22 +25,6 @@ a three-dimensional array containing normalized theta weights. Columns and rows This function estimates the weights of each theta (number of transmission events separating cases at two time points). A randomized transmission tree is drawn and the number of transmission events separating cases at two time points is calculated based on probabilies found in the Wallinga-Teunis matrix. } -\examples{ -case.times <- c(1,2,2,3,3) -gen <- c(0, 2/3, 1/3, 0, 0) -t.density <- gen/sum(gen) - -gen.time <- 2 # mean generation time - -wt <- est.wt.matrix(case.times=case.times, gen.t.dist=t.density) - -ngen <- round((max(case.times) - min(case.times)) / gen.time) + 1 # Number of generations - -a <- get.transdist.theta(wal.teun.mat=wt, - cases=case.times, - gen.t.mean=gen.time, - max.sep=ngen*2) -} \references{ Salje H, Cummings DAT and Lessler J (2016). “Estimating infectious disease transmission distances using the overall distribution of cases.” Epidemics, 17, pp. 10–18. ISSN 1755-4365, doi: \href{https://www.sciencedirect.com/science/article/pii/S1755436516300317}{10.1016/j.epidem.2016.10.001}. } @@ -54,3 +38,4 @@ Other transdist: \code{\link{est.transdist.bootstrap.ci}}, \author{ Justin Lessler, Henrik Salje, and John Giles } +\concept{transdist} diff --git a/man/sim.epidemic.Rd b/man/sim.epidemic.Rd index cd7fb8b..ab05e60 100644 --- a/man/sim.epidemic.Rd +++ b/man/sim.epidemic.Rd @@ -4,8 +4,8 @@ \alias{sim.epidemic} \title{Simulation of an epidemic in space and time} \usage{ -sim.epidemic(R, gen.t.mean, gen.t.sd, trans.kern.func, tot.generations = 10, - min.cases = 0, max.try = 1000) +sim.epidemic(R, gen.t.mean, gen.t.sd, trans.kern.func, + tot.generations = 10, min.cases = 0, max.try = 1000) } \arguments{ \item{R}{a scalar or a vector of length \code{tot.generations} providing the reproductive number for the epidemic. @@ -33,36 +33,6 @@ a function describing the spatial transmission kernel (\code{trans.kern.func}), of the generation time distribution (\code{gen.t.mean} and \code{gen.t.sd}) for the infecting pathogen. The function returns the location (\code{x}, \code{y}) and time (\code{t}) for each case of infection in the simulation. } -\examples{ -set.seed(1) - -dist.func <- alist(n=1, a=1/100, rexp(n, a)) # Exponential transmission kernel with mean = sd = 100 - -# Simulate epidemic with constant R value -a <- sim.epidemic(R=1.5, - gen.t.mean=7, - gen.t.sd=2, - tot.generations=15, - min.cases=100, - trans.kern.func=dist.func) - -sim.plot(a) - -# Simulate an epidemic with variable R value -r1 <- 2 -r2 <- 0.25 -tg <- 25 -R.vec <- seq(r1, r2, (r2 -r1)/(tg - 1)) - -b <- sim.epidemic(R=R.vec, - gen.t.mean=7, - gen.t.sd=2, - tot.generations=tg, - min.cases=100, - trans.kern.func=dist.func) - -sim.plot(b) -} \author{ Justin Lessler, Henrik Salje, and John Giles } diff --git a/src/spatialfuncs.c b/src/spatialfuncs.c index f952232..2c09e30 100644 --- a/src/spatialfuncs.c +++ b/src/spatialfuncs.c @@ -1,520 +1,520 @@ -#include -#include -#include -#include -#include - - -/****************************************************/ -/* convienence method for evaluating an R function */ -/* that returns 1,2 or 3 and takes in two rows of a */ -/* matrix. */ -/* @param Rfun the function */ -/* @param Rvect1 the first vector for the comparison*/ -/* @param Rvect2 the second vector for th comparison*/ -/**/ -/* @return the return value of Rfun */ -/****************************************************/ -double run_fun(SEXP Rfun, SEXP Rvect1, SEXP Rvect2) { - SEXP e, result; - - PROTECT(e = lang3(Rfun, Rvect1, Rvect2)); - result = eval(e, R_GlobalEnv); - UNPROTECT(1); - - return (REAL(result)[0]); -} - - -/*******************************************************************/ -/* convienence method to extract a row from a matrix as an R vector*/ -/* */ -/* @param matrix the matrix to extract from */ -/* @param row the row to extract */ -/* @return an SEXP pointing to an R vector */ -/*******************************************************************/ -SEXP extract_row(SEXP matrix, int row) { - int i; - SEXP rc; - SEXP dim = getAttrib(matrix, R_DimSymbol); - int rows = INTEGER(dim)[0]; - int cols = INTEGER(dim)[1]; - - PROTECT(rc = allocVector(REALSXP, cols)); - for (i=0; ir[i]) | (distr[i]) | (dist=r_low[i])) denom_cnt++; - - if (type[k] != *typeB) continue; - if ((dist<=r[i]) & (dist>=r_low[i])) num_cnt++; - } - } - - } else { - Rprintf("To be implemented\n"); - return; - } - //Rprintf("%d/%d\n",num_cnt,denom_cnt);//DEBUG - rc[i] = (double)num_cnt/denom_cnt; - } -} - - -/********************************************************************/ -/*theta function optimized for iterating through a list of two types*/ -/*finding ones that fulfill the distance requirement. */ -/* */ -/*@param type the type column */ -/*@param x the x coordinate */ -/*@param y the y coordinate */ -/*@param len the length of the three data arrays */ -/*@param typeA the "from" type */ -/*@param typeB the "to" type */ -/*@param r_low the low end of the range of values to look at */ -/*@param r the sequence of upper distances to consider */ -/*@param len_r the number of different Rs to consider */ -/* usually will be just the indicies */ -/*@param inds the indices into the original array, helps boot */ -/*@param rc the array of values that we are going to return */ -/********************************************************************/ -void get_theta_typed (int *type, - double *x, - double *y, - int *len, - int *typeA, - int *typeB, - double *r_low, - double *r, - int *len_r, - int *inds, - double *rc) { - - int i,j,k; - long long num_cnt, denom_cnt; /*counters for those filling conditions*/ - double dist; - - /*repeat the calculation for all r*/ - for (i=0;i<*len_r;i++) { - //Rprintf("%1.1f,%1.1f\n", r_low[i],r[i]); //DEBUG - /*zero out counts*/ - num_cnt = 0; - denom_cnt = 0; - - if (*typeA != -1) { - - for (j=0;j<*len;j++) { - if (type[j] != *typeA) continue; - - for (k=0;k<*len;k++) { - /*ignore pairs of the same observation*/ - if (inds[k]==inds[j]) continue; - - dist = sqrt(pow(x[j]-x[k],2)+pow(y[j]-y[k],2)); - - if ((dist<=r[i]) & (dist>=r_low[i])) { - if (type[k] == *typeB) { - num_cnt++; - } else { - denom_cnt++; - } - } - } - } - - } else { - Rprintf("To be implemented\n"); - return; - } - //Rprintf("%d/%d\n",num_cnt,denom_cnt);//DEBUG - rc[i] = (double)num_cnt/denom_cnt; - } -} - - - - - -/***********************************************************************/ -/*tau function for generic functions */ -/* */ -/* @param Rpostmat the matrix with the data in it. */ -/* @param Rfun the function to evaluate the relation between points */ -/* @param Rr the maximum distances to look at */ -/* @param Rr_min the minimum distances */ -/* @param Rcomparison_type 0 for represnetative, 1 for unrelated */ -/* @param Rinds indices into the original array, to help with bootstrapping*/ -/* @param Rxcol the column containing the x coordinate */ -/* @param Rycol the column containing the y coordinate */ -/***********************************************************************/ -SEXP get_tau (SEXP Rpostmat, - SEXP Rfun, - SEXP Rr, - SEXP Rr_low, - SEXP Rcomparison_type, - SEXP Rinds, - SEXP Rxcol, - SEXP Rycol) { - - int i; - SEXP divisor; - SEXP rc; - SEXP highR; - SEXP lowR; - int comparison_type = asInteger(Rcomparison_type); - - - PROTECT(highR=allocVector(REALSXP, 1)); - PROTECT(lowR=allocVector(REALSXP, 1)); - REAL(highR)[0] = R_PosInf; - REAL(lowR)[0] = 0; - - if (comparison_type==0) { - PROTECT(divisor=get_pi(Rpostmat,Rfun, highR, lowR, Rinds, Rxcol, Rycol)); - PROTECT(rc=get_pi(Rpostmat, Rfun, Rr, Rr_low, Rinds, Rxcol, Rycol)); - } else if (comparison_type==1) { - PROTECT(divisor=get_theta(Rpostmat,Rfun, highR, lowR, Rinds, Rxcol, Rycol)); - PROTECT(rc=get_theta(Rpostmat, Rfun, Rr, Rr_low, Rinds, Rxcol, Rycol)); - } else { - error("Invalid comparison_type."); - } - - for (i=0; i +#include +#include +#include +#include + + +/****************************************************/ +/* convienence method for evaluating an R function */ +/* that returns 1,2 or 3 and takes in two rows of a */ +/* matrix. */ +/* @param Rfun the function */ +/* @param Rvect1 the first vector for the comparison*/ +/* @param Rvect2 the second vector for th comparison*/ +/**/ +/* @return the return value of Rfun */ +/****************************************************/ +double run_fun(SEXP Rfun, SEXP Rvect1, SEXP Rvect2) { + SEXP e, result; + + PROTECT(e = lang3(Rfun, Rvect1, Rvect2)); + result = eval(e, R_GlobalEnv); + UNPROTECT(1); + + return (REAL(result)[0]); +} + + +/*******************************************************************/ +/* convienence method to extract a row from a matrix as an R vector*/ +/* */ +/* @param matrix the matrix to extract from */ +/* @param row the row to extract */ +/* @return an SEXP pointing to an R vector */ +/*******************************************************************/ +SEXP extract_row(SEXP matrix, int row) { + int i; + SEXP rc; + SEXP dim = getAttrib(matrix, R_DimSymbol); + int rows = INTEGER(dim)[0]; + int cols = INTEGER(dim)[1]; + + PROTECT(rc = allocVector(REALSXP, cols)); + for (i=0; i=r[i]) | (dist=r[i]) | (dist=r_low[i])) denom_cnt++; + + if (type[k] != *typeB) continue; + if ((dist=r_low[i])) num_cnt++; + } + } + + } else { + Rprintf("To be implemented\n"); + return; + } + //Rprintf("%d/%d\n",num_cnt,denom_cnt);//DEBUG + rc[i] = (double)num_cnt/denom_cnt; + } +} + + +/********************************************************************/ +/*theta function optimized for iterating through a list of two types*/ +/*finding ones that fulfill the distance requirement. */ +/* */ +/*@param type the type column */ +/*@param x the x coordinate */ +/*@param y the y coordinate */ +/*@param len the length of the three data arrays */ +/*@param typeA the "from" type */ +/*@param typeB the "to" type */ +/*@param r_low the low end of the range of values to look at */ +/*@param r the sequence of upper distances to consider */ +/*@param len_r the number of different Rs to consider */ +/* usually will be just the indicies */ +/*@param inds the indices into the original array, helps boot */ +/*@param rc the array of values that we are going to return */ +/********************************************************************/ +void get_theta_typed (int *type, + double *x, + double *y, + int *len, + int *typeA, + int *typeB, + double *r_low, + double *r, + int *len_r, + int *inds, + double *rc) { + + int i,j,k; + long long num_cnt, denom_cnt; /*counters for those filling conditions*/ + double dist; + + /*repeat the calculation for all r*/ + for (i=0;i<*len_r;i++) { + //Rprintf("%1.1f,%1.1f\n", r_low[i],r[i]); //DEBUG + /*zero out counts*/ + num_cnt = 0; + denom_cnt = 0; + + if (*typeA != -1) { + + for (j=0;j<*len;j++) { + if (type[j] != *typeA) continue; + + for (k=0;k<*len;k++) { + /*ignore pairs of the same observation*/ + if (inds[k]==inds[j]) continue; + + dist = sqrt(pow(x[j]-x[k],2)+pow(y[j]-y[k],2)); + + if ((dist=r_low[i])) { + if (type[k] == *typeB) { + num_cnt++; + } else { + denom_cnt++; + } + } + } + } + + } else { + Rprintf("To be implemented\n"); + return; + } + //Rprintf("%d/%d\n",num_cnt,denom_cnt);//DEBUG + rc[i] = (double)num_cnt/denom_cnt; + } +} + + + + + +/***********************************************************************/ +/*tau function for generic functions */ +/* */ +/* @param Rpostmat the matrix with the data in it. */ +/* @param Rfun the function to evaluate the relation between points */ +/* @param Rr the maximum distances to look at */ +/* @param Rr_min the minimum distances */ +/* @param Rcomparison_type 0 for represnetative, 1 for unrelated */ +/* @param Rinds indices into the original array, to help with bootstrapping*/ +/* @param Rxcol the column containing the x coordinate */ +/* @param Rycol the column containing the y coordinate */ +/***********************************************************************/ +SEXP get_tau (SEXP Rpostmat, + SEXP Rfun, + SEXP Rr, + SEXP Rr_low, + SEXP Rcomparison_type, + SEXP Rinds, + SEXP Rxcol, + SEXP Rycol) { + + int i; + SEXP divisor; + SEXP rc; + SEXP highR; + SEXP lowR; + int comparison_type = asInteger(Rcomparison_type); + + + PROTECT(highR=allocVector(REALSXP, 1)); + PROTECT(lowR=allocVector(REALSXP, 1)); + REAL(highR)[0] = R_PosInf; + REAL(lowR)[0] = 0; + + if (comparison_type==0) { + PROTECT(divisor=get_pi(Rpostmat,Rfun, highR, lowR, Rinds, Rxcol, Rycol)); + PROTECT(rc=get_pi(Rpostmat, Rfun, Rr, Rr_low, Rinds, Rxcol, Rycol)); + } else if (comparison_type==1) { + PROTECT(divisor=get_theta(Rpostmat,Rfun, highR, lowR, Rinds, Rxcol, Rycol)); + PROTECT(rc=get_theta(Rpostmat, Rfun, Rr, Rr_low, Rinds, Rxcol, Rycol)); + } else { + error("Invalid comparison_type."); + } + + for (i=0; i