Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@ S3method(print,ggwrm)
S3method(print,gtdrm)
S3method(print,gwavgm)
S3method(print,gwcorrm)
S3method(print,gwdam)
S3method(print,gwlcrm)
S3method(print,gwpcam)
S3method(print,gwrm)
Expand All @@ -39,6 +40,7 @@ export(gtdr_config)
export(gwaverage)
export(gwcorr_config)
export(gwcorrelation)
export(gwda)
export(gwpca)
export(gwr_basic)
export(gwr_generalized)
Expand Down
4 changes: 4 additions & 0 deletions R/RcppExports.R
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,10 @@ gw_correlation_cal <- function(x1, x2, coords, bw, adaptive, kernel, longlat, p,
.Call(`_GWmodel3_gw_correlation_cal`, x1, x2, coords, bw, adaptive, kernel, longlat, p, theta, initial_type, optim_bw_criterion, parallel_type, parallel_arg, variable_names, verbose)
}

gwda_cal <- function(x, y, coords, bw, adaptive, kernel, longlat, p, theta, method, parallel_type, parallel_arg) {
.Call(`_GWmodel3_gwda_cal`, x, y, coords, bw, adaptive, kernel, longlat, p, theta, method, parallel_type, parallel_arg)
}

gwpca_cal <- function(x, coords, bw, adaptive, kernel, longlat, p, theta, keep_components) {
.Call(`_GWmodel3_gwpca_cal`, x, coords, bw, adaptive, kernel, longlat, p, theta, keep_components)
}
Expand Down
180 changes: 180 additions & 0 deletions R/gwda.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,180 @@
#' Geographically Weighted Discriminant Analysis
#'
#' @param formula Regresison model.
#' @param data A `sf` objects.
#' @param bw Bandwidth value
#' @param adaptive Whether the bandwidth value is adaptive or not.
#' @param quantile Whether to calculate local quantiles.
#' @param kernel Kernel function used.
#' @param longlat Whether the coordinates
#' @param p Power of the Minkowski distance,
#' default to 2, i.e., Euclidean distance.
#' @param theta Angle in radian to roate the coordinate system, default to 0.
#' @param method Discriminant method for GWDA, can be set as
#' - `wqda`, weighted quadratic discriminant analysis
#' - `wlda`, weighted linear discriminant analysis
#' If not provided, this parameter will be set to `wqda`.
#' @param parallel_method Parallel method.
#' @param parallel_arg Parallel method argument.
#'
#' @return A `gwdam` object.
#'
#' @export
gwda <- function(
formula,
data,
bw,
adaptive = FALSE,
kernel = c("gaussian", "exp", "bisquare", "tricube", "boxcar"),
longlat = FALSE,
p = 2.0,
theta = 0.0,
method = c("wqda", "wlda"),
parallel_method = c("no", "omp"),
parallel_arg = c(0))
{
### Check args
kernel <- match.arg(kernel)
parallel_method <- match.arg(parallel_method)
attr(data, "na.action") <- getOption("na.action")

### Extract coords
data <- do.call(na.action(data), args = list(data))
coords <- as.matrix(sf::st_coordinates(sf::st_centroid(data)))
if (is.null(coords) || nrow(coords) != nrow(data)) {
stop("Missing coordinates.")
}

### Extract variables
mc <- match.call(expand.dots = FALSE)
mf <- model.frame(formula = formula(update(formula, ~ . - 1)), data = sf::st_drop_geometry(data))
mt <- attr(mf, "terms")
y <- model.extract(mf, "response")
y <- as.character(y)
x <- model.matrix(mt, mf)
dep_var <- as.character(attr(terms(mf), "variables")[[2]])
indep_vars <- colnames(x)

### Check whether bandwidth is valid.
if (missing(bw)) {
stop("Bandwidth is missing.")
} else if (is.numeric(bw) || is.integer(bw)) {
if (bw == 0) stop("Bandwidth must be non-zero.")
} else {
stop("Bandwidth must be a number.")
}

method <- match.arg(method)
wqda_method <- ifelse(method == "wqda", TRUE, FALSE)

c_results <- gwda_cal(
x, y, coords, bw,
adaptive, enum(kernel),
longlat, p, theta,
wqda_method,
enum_list(parallel_method, parallel_types), parallel_arg
)

sdf_data <- data.frame()
res <- c_results$res
group <- as.data.frame(c_results$group)
probs <- c_results$probs
pmax <- c_results$pmax
entropy <- c_results$entropy

colnames(res) <- indep_vars
colnames(group) <- "groups"
colnames(probs) <- paste0(indep_vars, "_probs")
colnames(pmax) <- "pmax"
colnames(entropy) <- "entropy"
sdf_data <- as.data.frame(cbind(res, group, probs, pmax, entropy))

sdf_data$geometry <- sf::st_geometry(data)
sdf <- sf::st_sf(sdf_data)

### Return result
gwdam <- list(
SDF = sdf,
args = list(
x = x,
y = y,
correctRate = c_results$correctRate,
coords = coords,
bw = bw,
adaptive = adaptive,
kernel = kernel,
longlat = longlat,
p = p,
theta = theta,
method = method,
parallel_method = parallel_method,
parallel_arg = parallel_arg
),
call = mc,
indep_vars = indep_vars,
dep_var = dep_var
)
class(gwdam) <- "gwdam"
gwdam
}

#' Print description of a `gwdam` object
#'
#' @param x An `gwdam` object returned by [gwr_basic()].
#' @param decimal_fmt The format string passing to [base::sprintf()].
#' @inheritDotParams print_table_md
#'
#' @method print gwdam
#' @rdname print
#' @export
print.gwdam <- function(x, ..., decimal_fmt) {
if (!inherits(x, "gwdam")) {
stop("It's not a 'gwdam' object.")
}

cat(" ***********************************************************************\n")
cat(" * Package GWmodel3 *\n")
cat(" ***********************************************************************\n")
cat(" * Results of Geographically Discriminant Analysis *\n")
cat(" ***********************************************************************\n")
cat("\n *********************Model Calibration Information*********************\n")

cat(" Formula:", deparse(x$call$formula), fill = T)
cat(" Data:", deparse(x$call$data), fill = T)
cat(" Method:", x$args$method, fill = T)
cat(" Number of summary points:", nrow(x$args$x), fill = T)
cat(" Kernel function:", x$args$kernel, fill = T)
cat(" Bandwidth:", x$args$bw,
ifelse(x$args$adaptive, "(Nearest Neighbours)", "(Meters)"),
fill = T
)
distance_type <- "Euclidean"
if (x$args$longlat) {
distance_type <- "Geodetic"
} else if (x$args$p == 2) {
distance_type <- "Euclidean"
} else if (x$args$p == 1) {
distance_type <- "Manhattan"
} else if (is.infinite(x$args$p)) {
distance_type <- "Chebyshev"
} else {
distance_type <- "Generalized Minkowski"
}
distance_rotated <- (x$args$theta != 0 && x$args$p != 2 && !x$args$longlat)
cat(" Distance:", distance_type, ifelse(distance_rotated, " (rotated)", ""), fill = T)
res <- st_drop_geometry(x$SDF)
group_counts <- as.data.frame(table(res$groups))
cat(" Discrimination Result:", fill = T)
colnames(group_counts) <- c("Group", "Count")
output <- capture.output(print(group_counts, row.names = FALSE))
output <- paste0(" ", output)
cat(output, sep = "\n")

cat("\n ***********************Local Summary Statistics************************", fill = T)
numeric_res <- res[sapply(res, is.numeric)]
df <- as.data.frame(t(sapply(numeric_res, summary)))
output <- capture.output(print(df))
output <- paste0(" ", output)
cat(output, sep = "\n")
cat(" ***********************************************************************\n")
}
57 changes: 57 additions & 0 deletions man/gwda.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

11 changes: 8 additions & 3 deletions man/print.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

2 changes: 2 additions & 0 deletions src/Makevars.in
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,7 @@ OBJECTS_LIBGWMODEL = \
libgwmodel/src/gwmodelpp/GWRGeneralized.o \
libgwmodel/src/gwmodelpp/GWRMultiscale.o \
libgwmodel/src/gwmodelpp/GWAverage.o \
libgwmodel/src/gwmodelpp/GWDA.o \
libgwmodel/src/gwmodelpp/GWCorrelation.o \
libgwmodel/src/gwmodelpp/GWRLocalCollinearity.o \
libgwmodel/src/gwmodelpp/GWPCA.o \
Expand All @@ -73,6 +74,7 @@ OBJECTS_GWMODEL = \
gwr_localcollinearity.o \
gtdr.o \
gwaverage.o \
gwda.o \
gwcorrelation.o \
gwpca.o \
RcppExports.o
Expand Down
2 changes: 2 additions & 0 deletions src/Makevars.win
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,7 @@ OBJECTS_LIBGWMODEL = \
libgwmodel/src/gwmodelpp/GWRGeneralized.o \
libgwmodel/src/gwmodelpp/GWRMultiscale.o \
libgwmodel/src/gwmodelpp/GWAverage.o \
libgwmodel/src/gwmodelpp/GWDA.o \
libgwmodel/src/gwmodelpp/GWCorrelation.o \
libgwmodel/src/gwmodelpp/GWRLocalCollinearity.o \
libgwmodel/src/gwmodelpp/GWPCA.o \
Expand All @@ -75,6 +76,7 @@ OBJECTS_GWMODEL = \
gwr_localcollinearity.o \
gtdr.o \
gwaverage.o \
gwda.o \
gwcorrelation.o \
gwpca.o \
RcppExports.o
Expand Down
23 changes: 23 additions & 0 deletions src/RcppExports.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -86,6 +86,28 @@ BEGIN_RCPP
return rcpp_result_gen;
END_RCPP
}
// gwda_cal
List gwda_cal(const arma::mat& x, std::vector<std::string>& y, const arma::mat& coords, double bw, bool adaptive, size_t kernel, bool longlat, double p, double theta, bool method, size_t parallel_type, const IntegerVector& parallel_arg);
RcppExport SEXP _GWmodel3_gwda_cal(SEXP xSEXP, SEXP ySEXP, SEXP coordsSEXP, SEXP bwSEXP, SEXP adaptiveSEXP, SEXP kernelSEXP, SEXP longlatSEXP, SEXP pSEXP, SEXP thetaSEXP, SEXP methodSEXP, SEXP parallel_typeSEXP, SEXP parallel_argSEXP) {
BEGIN_RCPP
Rcpp::RObject rcpp_result_gen;
Rcpp::RNGScope rcpp_rngScope_gen;
Rcpp::traits::input_parameter< const arma::mat& >::type x(xSEXP);
Rcpp::traits::input_parameter< std::vector<std::string>& >::type y(ySEXP);
Rcpp::traits::input_parameter< const arma::mat& >::type coords(coordsSEXP);
Rcpp::traits::input_parameter< double >::type bw(bwSEXP);
Rcpp::traits::input_parameter< bool >::type adaptive(adaptiveSEXP);
Rcpp::traits::input_parameter< size_t >::type kernel(kernelSEXP);
Rcpp::traits::input_parameter< bool >::type longlat(longlatSEXP);
Rcpp::traits::input_parameter< double >::type p(pSEXP);
Rcpp::traits::input_parameter< double >::type theta(thetaSEXP);
Rcpp::traits::input_parameter< bool >::type method(methodSEXP);
Rcpp::traits::input_parameter< size_t >::type parallel_type(parallel_typeSEXP);
Rcpp::traits::input_parameter< const IntegerVector& >::type parallel_arg(parallel_argSEXP);
rcpp_result_gen = Rcpp::wrap(gwda_cal(x, y, coords, bw, adaptive, kernel, longlat, p, theta, method, parallel_type, parallel_arg));
return rcpp_result_gen;
END_RCPP
}
// gwpca_cal
List gwpca_cal(const arma::mat& x, const arma::mat& coords, double bw, bool adaptive, size_t kernel, bool longlat, double p, double theta, size_t keep_components);
RcppExport SEXP _GWmodel3_gwpca_cal(SEXP xSEXP, SEXP coordsSEXP, SEXP bwSEXP, SEXP adaptiveSEXP, SEXP kernelSEXP, SEXP longlatSEXP, SEXP pSEXP, SEXP thetaSEXP, SEXP keep_componentsSEXP) {
Expand Down Expand Up @@ -279,6 +301,7 @@ static const R_CallMethodDef CallEntries[] = {
{"_GWmodel3_gtdr_fit", (DL_FUNC) &_GWmodel3_gtdr_fit, 19},
{"_GWmodel3_gwaverage_fit", (DL_FUNC) &_GWmodel3_gwaverage_fit, 11},
{"_GWmodel3_gw_correlation_cal", (DL_FUNC) &_GWmodel3_gw_correlation_cal, 15},
{"_GWmodel3_gwda_cal", (DL_FUNC) &_GWmodel3_gwda_cal, 12},
{"_GWmodel3_gwpca_cal", (DL_FUNC) &_GWmodel3_gwpca_cal, 9},
{"_GWmodel3_gwr_basic_fit", (DL_FUNC) &_GWmodel3_gwr_basic_fit, 22},
{"_GWmodel3_gwr_basic_predict", (DL_FUNC) &_GWmodel3_gwr_basic_predict, 14},
Expand Down
Loading
Loading