[spcopula-commits] r102 - / pkg pkg/R pkg/man pkg/src
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Sep 1 13:10:27 CEST 2013
Author: mappe
Date: 2013-09-01 13:10:27 +0200 (Sun, 01 Sep 2013)
New Revision: 102
Added:
pkg/man/getNeighbours.experimental.Rd
pkg/src/
pkg/src/Rinterface.h
pkg/src/spatialPreparation.c
thoughts-Marius/
Modified:
pkg/NAMESPACE
pkg/R/spatialPreparation.R
Log:
- added thoughts-Marius folder
- added C implementation of getNeighbours() in getNeighbours.experimental()
Modified: pkg/NAMESPACE
===================================================================
--- pkg/NAMESPACE 2013-08-01 10:31:15 UTC (rev 101)
+++ pkg/NAMESPACE 2013-09-01 11:10:27 UTC (rev 102)
@@ -27,7 +27,7 @@
# export(setSizeLim)
# spatial
-export(getNeighbours,getStNeighbours)
+export(getNeighbours,getNeighbours.experimental, getStNeighbours)
export(calcBins)
export(calcSpTreeDists, dropSpTree)
Modified: pkg/R/spatialPreparation.R
===================================================================
--- pkg/R/spatialPreparation.R 2013-08-01 10:31:15 UTC (rev 101)
+++ pkg/R/spatialPreparation.R 2013-09-01 11:10:27 UTC (rev 102)
@@ -98,6 +98,60 @@
predLocs, prediction, var))
}
+
+
+getNeighbours.experimental = function (dataLocs, predLocs, var = names(dataLocs)[1], size = 5,
+ prediction = FALSE, min.dist = 0.01)
+{
+ stopifnot((!prediction && missing(predLocs)) || (prediction &&
+ !missing(predLocs)))
+ stopifnot(min.dist > 0 || prediction)
+ if (missing(predLocs) && !prediction)
+ predLocs = dataLocs
+ stopifnot(is(predLocs, "Spatial"))
+ if (any(is.na(match(var, names(dataLocs)))))
+ stop("At least one of the variables is unkown or is not part of the data.")
+ nLocs <- length(predLocs)
+ size <- min(size, length(dataLocs) + prediction)
+
+ allLocs <- matrix(0, nLocs, size)
+ allDists <- matrix(0, nLocs, size - 1)
+ result = .C("getNeighbours_2d",
+ allLocs = as.double(allLocs),
+ allDists = as.double(allDists),
+ as.double(coordinates(dataLocs)),
+ as.double(coordinates(predLocs)),
+ as.integer(length(dataLocs)),
+ as.integer(length(predLocs)),
+ as.double(min.dist),
+ as.integer(size),
+ as.integer(prediction),
+ PACKAGE="spcopula")
+
+ if (!prediction) allData <- matrix(dataLocs[result$allLocs, var, drop = F]@data[[1]], nLocs, size)
+ else allData <- matrix(c(rep(NA,nLocs),dataLocs[result$allLocs[(nLocs+1):(nLocs*size)], var, drop = F]@data[[1]]), nLocs, size)
+
+ if (!prediction)
+ predLocs <- NULL
+
+ colnames(allData) <- paste(paste("N", rep(0:(size - 1), each = length(var)),
+ sep = ""), rep(var, size), sep = ".")
+
+
+ return(neighbourhood(allData, matrix(result$allDists,nLocs, size - 1), matrix(result$allLocs, nLocs, size), dataLocs,
+ predLocs, prediction, var))
+}
+
+
+
+
+
+
+
+
+
+
+
#############
## BINNING ##
#############
Added: pkg/man/getNeighbours.experimental.Rd
===================================================================
--- pkg/man/getNeighbours.experimental.Rd (rev 0)
+++ pkg/man/getNeighbours.experimental.Rd 2013-09-01 11:10:27 UTC (rev 102)
@@ -0,0 +1,47 @@
+\name{getNeighbours.experimental}
+\alias{getNeighbours.experimental}
+
+\title{
+Creating Local Neighbourhoods
+}
+\description{
+This function calculates a local neighbourhood to be used for fitting of spatial/spatio-temporal vine copulas and for prediction using spatial/spatio-temporal vine copulas.
+}
+\usage{
+getNeighbours.experimental(dataLocs, predLocs, var = names(dataLocs)[1], size = 5, prediction=FALSE, min.dist = 0.01)
+}
+\arguments{
+ \item{dataLocs}{
+some spatial data frame holding the data used for estimation/prediction
+}
+ \item{predLocs}{
+A spatial object defining the prediction locations, might be missing if the neighbourhood is used for fitting.
+}
+ \item{var}{
+the variable name of interest, by default the first variable is used.
+}
+ \item{size}{
+The size of the neighbourhood including the location of interest (for fitting as well for prediction).
+}
+ \item{prediction}{whether the neighbourhood should be used for prediction (TRUE) or spatial/Spatio-temporal vine copula fitting.}
+ \item{min.dist}{
+the minimal distance for a location to be included. Must be larger than 0 for fitting purposes and might be 0 for prediction.
+}
+}
+\value{
+An object of \code{\linkS4class{neighbourhood}}.
+}
+\author{
+Benedikt Graeler
+}
+
+\seealso{
+See \code{\link{neighbourhood}} for the native constructor of a \code{\linkS4class{neighbourhood}} class.
+}
+\examples{
+spdf <- data.frame(x=c(112,154,212,289,345),y=c(124,198,85,168,346),measure=rlnorm(5))
+coordinates(spdf) <- ~x+y
+
+getNeighbours.experimental(spdf,size=4)
+}
+\keyword{ spatial }
\ No newline at end of file
Added: pkg/src/Rinterface.h
===================================================================
--- pkg/src/Rinterface.h (rev 0)
+++ pkg/src/Rinterface.h 2013-09-01 11:10:27 UTC (rev 102)
@@ -0,0 +1,33 @@
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+
+
+#ifdef _WIN32
+#define R_EXPORT __declspec(dllexport)
+#else
+#define R_EXPORT
+#endif
+
+
+
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+/* FORWARD DECLARATIONS */
+
+void R_EXPORT distMatrix(double *out, double *xa, double *xb, int *na, int*nb, int *d);
+void R_EXPORT distMatrix_2d(double *out, double *xa, double *xb, int *na, int*nb);
+void R_EXPORT distMatrixInner(double *out, double *xa, int *na, int *d);
+void R_EXPORT distMatrixInner_2d(double *out, double *xa, int *na);
+
+void R_EXPORT getNeighbours_2d(double *out_locs, double *out_dists, double *data_locs, double *pred_locs, int *n_data, int *n_pred, double *min_dist, int *k, int *prediction);
+
+
+#ifdef __cplusplus
+}
+#endif
\ No newline at end of file
Added: pkg/src/spatialPreparation.c
===================================================================
--- pkg/src/spatialPreparation.c (rev 0)
+++ pkg/src/spatialPreparation.c 2013-09-01 11:10:27 UTC (rev 102)
@@ -0,0 +1,244 @@
+
+#include "Rinterface.h"
+
+
+
+
+
+/* Helper Functions */
+inline double dist_2d(double ax, double ay, double bx, double by) {
+ return(sqrt( (ax-bx)*(ax-bx) + (ay-by)*(ay-by) ));
+}
+
+
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+
+
+
+
+
+/**
+* Computes the distance matrix of two sets of n-dimensional points
+* @param out result (na,nb) matrix (mast be allocated in R or in calling function before)
+* @param xa coordinates of the 1st set of points, each row of this matrix is a point
+* @param xb coordinates of the 2nd set of points, each row of this matrix is a point
+* @param na number of points in xa
+* @param nb number of points in xb
+* @param d number of dimensions, equals the number of columns in xa and xb
+**/
+void R_EXPORT distMatrix(double *out, double *xa, double *xb, int *na, int*nb, int *d) {
+ if (*d == 2) {
+ distMatrix_2d(out, xa, xb, na, nb);
+ }
+ else
+ {
+ for (int i=0;i<*na; ++i) {
+ for (int j=0;j<*nb; ++j) {
+ out[j * (*na) + i] = 0;
+
+ for (int k=0; k<*d; ++k) {
+ out[j * (*na) + i] += (xa[k * (*na) + i] - xb[k * (*nb) + j]) * (xa[k * (*na) + i] - xb[k * (*nb) + j]); // In R: col major
+
+ }
+ out[j * (*na) + i] = sqrt(out[j * (*na) + i]);
+ }
+ }
+ }
+}
+
+
+
+/**
+* Computes the distance matrix of two sets of 2-dimensional points
+* @param out result (na,nb) matrix (mast be allocated in R or in calling function before)
+* @param xa coordinates of the 1st set of points, each row of this matrix is a point
+* @param xb coordinates of the 2nd set of points, each row of this matrix is a point
+* @param na number of points in xa
+* @param nb number of points in xb
+**/
+void R_EXPORT distMatrix_2d(double *out, double *xa, double *xb, int *na, int*nb) {
+ for (int i=0;i<*na; ++i) {
+ for (int j=0;j<*nb; ++j) {
+ out[j * (*na) + i] = dist_2d(xa[i],xa[*na + i], xb[j], xb[*nb + j]);
+ }
+ }
+}
+
+
+
+/**
+* Computes the "inner" distance matrix of one set of n-dimensional points
+* @param out result (na,na) matrix (mast be allocated in R or in calling function before)
+* @param xa coordinates of the 1st set of points, each row of this matrix is a point
+* @param na number of points in xa
+* @param d number of dimensions, equals the number of columns in xa and xb
+**/
+void R_EXPORT distMatrixInner(double *out, double *xa, int *na, int *d) {
+ if (*d == 2) {
+ distMatrixInner_2d(out, xa, na);
+ }
+ else
+ {
+ for (int i=0;i<*na; ++i) {
+ for (int j=i+1;j<*na; ++j) {
+ out[j * (*na) + i] = 0;
+ for (int k=0; k<*d; ++k) {
+ out[j * (*na) + i] += (xa[k * (*na) + i] - xa[k * (*na) + j]) * (xa[k * (*na) + i] - xa[k * (*na) + j]); // In R: col major
+
+ }
+ out[j * (*na) + i] = sqrt(out[j * (*na) + i]);
+ out[i * (*na) + j] = out[j * (*na) + i];
+ }
+ }
+ for (int i=0;i<*na; ++i) {
+ out[i * (*na) + i] = 0.0;
+ }
+ }
+}
+
+
+
+/**
+* Computes the "inner" distance matrix of one set of n-dimensional points
+* @param out result (na,na) matrix (mast be allocated in R or in calling function before)
+* @param xa coordinates of the 1st set of points, each row of this matrix is a point
+* @param na number of points in xa
+**/
+void R_EXPORT distMatrixInner_2d(double *out, double *xa, int *na) {
+ for (int i=0;i<*na; ++i) {
+ for (int j=i+1;j<*na; ++j) {
+ out[j * (*na) + i] = dist_2d(xa[i],xa[*na + i], xa[j], xa[*na + j]);
+ out[i * (*na) + j] = out[j * (*na) + i];
+ }
+ }
+ // set diagonal to zero
+ for (int i=0;i<*na; ++i) {
+ out[i * (*na) + i] = 0.0;
+ }
+}
+
+
+
+
+
+
+/*for (i in 1:nLocs) {
+ tempDists <- spDists(dataLocs, predLocs[i, ])
+ tempDists[tempDists < min.dist] <- Inf
+ spLocs <- order(tempDists)[1:(size - 1)]
+
+ allLocs[i,] <- c(i, spLocs)
+ allDists[i,] <- tempDists[spLocs]
+ allData[i,(prediction+1):size] <- dataLocs[c(i[!prediction],spLocs),
+ var, drop = F]@data[[1]]
+ }
+ */
+
+/**
+*
+*
+*
+*
+* @param k size of the neighbourhood including the location of interest
+*/
+void R_EXPORT getNeighbours_2d(double *out_locs, double *out_dists, double *data_locs, double *pred_locs, int *n_data, int *n_pred, double *min_dist, int *k, int *prediction) {
+ double *C = (double*)malloc((*n_data) * (*n_pred) * sizeof(double));
+ double *knn_dists = (double*)malloc((*k-1) * sizeof(double));
+ int *knn_idxs = (int*)malloc((*k-1) * sizeof(int));
+
+ // Depending on whether pred_locs are the same as data_locs, use optimized "inner" distance matrix calculation
+ if (!prediction) {
+ // Compute inner distance matrix and store in C
+ distMatrixInner_2d(C, data_locs, n_data);
+ }
+ else {
+
+ // Compute distance matrix and store in C
+ distMatrix_2d(C, data_locs, pred_locs, n_data, n_pred);
+ }
+
+ // for each pred location
+ for (int i=0; i<*n_pred; ++i) {
+ knn_dists[0] = INFINITY;
+ knn_idxs[0] = -1;
+
+ for (int j=0; j<*n_data; ++j) {
+ if (C[i* (*n_data) + j] < *min_dist) continue;
+ int nbrank = 0;
+ while (C[i* (*n_data) + j] > knn_dists[nbrank] && nbrank < (*k-1)) ++nbrank;
+ if (nbrank < (*k-1)) {
+ for (int s=(*k-3); s >= nbrank; --s) {
+ knn_dists[s + 1] = knn_dists[s]; // shift array if new element in knn array
+ knn_idxs[s + 1] = knn_idxs[s];
+ }
+ knn_dists[nbrank] = C[i* (*n_data) + j];
+ knn_idxs[nbrank] = j+1; // Indices in R start with 1
+ }
+ }
+
+ out_locs[i] = i+1; // Indices in R start with 1
+ for (int j=1; j<(*k); ++j) {
+ // Assign c(i, knn_idxs[1:size-1]) as i-th row in out_locs
+ out_locs[j*(*n_pred) + i] = knn_idxs[j-1];
+ // Assign knn_dists[1:size-1] as i-th row in out_dists
+ out_dists[(j-1)*(*n_pred) + i] = knn_dists[j-1];
+ }
+
+ }
+
+ free(knn_idxs);
+ free(knn_dists);
+ free(C);
+}
+
+
+
+
+
+
+/*
+calcSpLagInd <- function(data, boundaries) {
+ lags <- vector("list",length(boundaries))
+
+ dists <- spDists(data)
+ nlocs <- length(data)
+
+ for (i in 1:(nlocs-1)) {
+ for (j in (i+1):nlocs) {
+ d <- dists[i,j]
+ for ( k in 1:length(boundaries)) {
+ if (d < boundaries[k]) {
+ lags[[k]] <- rbind(lags[[k]],c(i,j,d))
+ break()
+ }
+ }
+ }
+ }
+ return(lags)
+}
+*/
+
+/**
+*
+* @
+*
+*
+*
+* OpenMP?
+*/
+// SEXP calcSpLagInd(double *o double *idata, double *iboundaries, int *in_boundaries) {
+
+
+ // for
+
+// }
+
+
+
+#ifdef __cplusplus
+}
+#endif
\ No newline at end of file
More information about the spcopula-commits
mailing list