/
predict.stratafit.R
88 lines (75 loc) · 3.39 KB
/
predict.stratafit.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
#' Perform Finite Population Block Kriging
#'
#' Uses an object of class \code{stratafit} from the \code{\link{stratafit}()}
#' function to predict the response on the unsampled sites for separate strata.
#' The column of the data set that has the response should have
#' numeric values for the observed response on the sampled sites and
#' `NA` for any site that was not sampled. Note that there is no
#' \code{newdata} argument to \code{predict.stratafit()}: any point in
#' space for which a prediction is needed should be included in the
#' original data set in \code{\link{stratafit}()} with the response
#' variable as \code{NA}.
#'
#' @param object is an object generated from \code{\link{stratafit}()}
#' @param wtscol is the name of the column that contains the weights
#' for prediction. The default setting predicts the population total
#' @param conf_level, by default, 0.90, this is the desired
#' confidence level for a prediction interval
#' @param ... further arguments passed to or from other methods.
#' @return a list with \itemize{
#' \item the estimated population total
#' \item the estimated prediction variance
#' \item a data frame containing \enumerate{
#' \item x-coordinates
#' \item y-coordinates
#' \item density predictions
#' \item count predictions
#' \item site-by-site density prediction variances
#' \item site-by-site count prediction variances
#' \item indicator variable for whether or not the
#' each site was sampled
#' \item estimated mean for each site
#' \item area of each site
#' }
#' \item vector with estimated covariance parameters
#' \item the formula used to fit the model in \code{slmfit()}
#' }
#' @examples
#' data(exampledataset) ## load a toy data set
#' exampledataset$strata <- c(rep("A", 19), rep("B", 21))
#' strataobj <- stratafit(formula = counts ~ pred1 + pred2,
#' data = exampledataset, stratacol = "strata",
#' xcoordcol = 'xcoords', ycoordcol = 'ycoords', areacol = 'areavar')
#' predict(strataobj)
#' @importFrom stats predict
#' @export
predict.stratafit <- function(object, wtscol = NULL,
conf_level = 0.9, ...) {
predict_outs <- lapply(object, FUN = predict, wtscol = wtscol)
preds <- lapply(predict_outs, `[[`, 1)
vars <- lapply(predict_outs, `[[`, 2)
dfs <- lapply(predict_outs, `[[`, 3)
pred_out <- do.call("sum", preds)
ses <- lapply(vars, sqrt)
totalvar_out <- do.call("sum", vars)
se_out <- sqrt(totalvar_out)
df_out <- do.call("rbind", dfs)
strata_out <- list(preds, pred_out, vars, se_out, df_out)
total_out <- matrix(c(pred_out, se_out), nrow = 1)
pred_se_mat <- rbind(cbind(matrix(unlist(preds)),
matrix(unlist(ses))), total_out)
lbs <- matrix(pred_se_mat[ ,1] + 1 *
stats::qnorm((1 - conf_level) / 2) * pred_se_mat[ ,2])
ubs <- matrix(pred_se_mat[ ,1] - 1 *
stats::qnorm((1 - conf_level) / 2) * pred_se_mat[ ,2])
output_mat <- cbind(pred_se_mat, lbs, ubs)
rownames(output_mat) <- c(names(dfs), "Total")
colnames(output_mat) <- c("Prediction", "SE",
paste(conf_level * 100, "% LB", sep = ""),
paste(conf_level * 100,
"% UB", sep = ""))
output_obj <- list(output_mat, df_out)
names(output_obj) <- c("summary_info", "Pred_df")
class(output_obj) <- "predict.stratafit"
return(output_obj)
}