|
DI offnw.probCalculate off-network probabilities. DESCRIPTION: Takes a network object and then uses Krig.nw (which in turn uses the Fields function, Krig, to calculate the probability of exceeding a given value on any grid point off of the network. This function is called by plot.prob USAGE: offnw.prob(nw.obj=new.network, Y, K=85, GCV=T, deg.f=NA, Rho=NA, Sigma2=NA, plts="Probabilities", ...) REQUIRED ARGUMENTS: nw.obj Network object created by make.network.obj Y Vector of measurements recorded at each point on the network. K Exceedence value; P[ Z > K]. Default is 85. OPTIONAL ARGUMENTS: GCV Use Generalized Cross-Validation criterion for kriging? deg.f degrees of freedom; for use with the Krig function. Rho correlation parameter for use with the Krig function. Sigma2 variance parameter for use with the Krig function. plts Which plot you wish to view. Choices are: probabilities, the means as found by kriging, the standard deviations as found by kriging, diagnostics of the Krig function or all of the above. This is more important for the function plot.prob which actually calls offnw.prob VALUE: Returns a list object giving the means, standard deviations, probabilities, etc... SEE ALSO: Krig, Krig.nw, plot.prob EXAMPLES: # The function is currently defined as function(nw.obj = new.network, Y, K = 85, GCV = T, deg.f = NA, Rho = NA, Sigma2 = NA, plts = "Probabilities", ...){ x0 <- make.surface.grid(nw.obj$grid) y0 <- Y if(is.null(nw.obj$cov.obj$mean.obj)) { mean.obj <- 0 class(mean.obj) <- "constant" } else if(is.numeric(nw.obj$cov.obj$mean.obj)) { mean.obj <- nw.obj$cov.obj$mean.obj class(mean.obj) <- "constant" } else mean.obj <- get(nw.obj$cov.obj$mean.obj) if(is.null(nw.obj$cov.obj$sd.obj)) { sd.obj <- 1 class(sd.obj) <- "constant" } else if(is.numeric(nw.obj$cov.obj$sd.obj)) { sd.obj <- nw.obj$cov.obj$sd.obj class(sd.obj) <- "constant" } else sd.obj <- get(nw.obj$cov.obj$sd.obj) if(GCV) Krig.obj <- Krig.nw(nw.obj, Y = y0, rho = Rho, sigma2 = Sigma2, meanobj = mean.obj, sdobj = sd.obj, ...) else Krig.obj <- Krig.nw(nw.obj, Y = y0, df = deg.f, rho = Rho, sigma2 = Sigma2, meanobj = mean.obj, sdobj = sd.obj, ...) mtemp <- predict(Krig.obj, x0) setemp <- predict.se(Krig.obj, x0) p <- 1 - pnorm((K - mtemp)/setemp) if(plts == "diagnostics" | plts == "all") { list(prob = as.surface(x0, p), mean = as.surface(x0, mtemp), se = as.surface(x0, setemp), K = K, fit = Krig.obj) } else { list(prob = as.surface(x0, p), mean = as.surface(x0, mtemp), se = as.surface(x0, setemp), K = K) } } |
This is software for statistical research and should not be used for commercial purposes. The authors do not guarantee the correctness of any function or program in this package.