diff --git a/.Rbuildignore b/.Rbuildignore
index 4b349a88..03475a99 100644
--- a/.Rbuildignore
+++ b/.Rbuildignore
@@ -13,3 +13,5 @@
 ^cran-comments\.md$
 ^pull_request_template$
 PULL_REQUEST_TEMPLATE.md
+.claude
+.idea
diff --git a/.gitignore b/.gitignore
index 573f549b..f4edd3a2 100644
--- a/.gitignore
+++ b/.gitignore
@@ -6,3 +6,5 @@
 .Rproj.user
 .DS_Store
 .Rapp.history
+.claude
+.idea
diff --git a/R/corDS.R b/R/corDS.R
index abc73145..6133b4e0 100644
--- a/R/corDS.R
+++ b/R/corDS.R
@@ -14,8 +14,10 @@
 #' sum of squares of each variable. The first disclosure control checks that the number of variables is
 #' not bigger than a percentage of the individual-level records (the allowed percentage is pre-specified
 #' by the 'nfilter.glm'). The second disclosure control checks that none of them is dichotomous with a
-#' level having fewer counts than the pre-specified 'nfilter.tab' threshold.
+#' level having fewer counts than the pre-specified 'nfilter.tab' threshold. The list also includes
+#' \code{class}, the class of the input object for client-side consistency checking.
 #' @author Paul Burton, and Demetris Avraam for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #'
 corDS <- function(x=NULL, y=NULL){
@@ -27,14 +29,21 @@ corDS <- function(x=NULL, y=NULL){
   nfilter.glm <- as.numeric(thr$nfilter.glm)
   #############################################################
   
-  x.val <- eval(parse(text=x), envir = parent.frame())
+  x.val <- .loadServersideObject(x)
+  .checkClass(obj = x.val, obj_name = x, permitted_classes = c("numeric", "integer", "matrix", "data.frame"))
+
   if (!is.null(y)){
-    y.val <- eval(parse(text=y), envir = parent.frame())
+    y.val <- .loadServersideObject(y)
+    .checkClass(obj = y.val, obj_name = y, permitted_classes = c("numeric", "integer", "matrix", "data.frame"))
   }
   else{
     y.val <- NULL
   }
-  
+
+  if (is.null(y.val) && any(class(x.val) %in% c("numeric", "integer"))) {
+    stop("If x is a numeric vector, y must also be a numeric vector.", call. = FALSE)
+  }
+
   # create a data frame for the variables
   if (is.null(y.val)){
     dataframe <- as.data.frame(x.val)
@@ -165,7 +174,7 @@ corDS <- function(x=NULL, y=NULL){
     
   }
   
-  return(list(sums.of.products=sums.of.products, sums=sums, complete.counts=complete.counts, na.counts=na.counts, sums.of.squares=sums.of.squares))
+  return(list(sums.of.products=sums.of.products, sums=sums, complete.counts=complete.counts, na.counts=na.counts, sums.of.squares=sums.of.squares, class=class(x.val)))
   
 }
 # AGGREGATE FUNCTION
diff --git a/R/corTestDS.R b/R/corTestDS.R
index ef5aac33..777bdce8 100644
--- a/R/corTestDS.R
+++ b/R/corTestDS.R
@@ -11,14 +11,18 @@
 #' @param conf.level confidence level for the returned confidence interval. Currently
 #' only used for the Pearson product moment correlation coefficient if there are at least
 #' 4 complete pairs of observations.
-#' @return the results of the correlation test.
+#' @return a list with the results of the correlation test and \code{class}, the class of the
+#' input object for client-side consistency checking.
 #' @author Demetris Avraam, for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #'
 corTestDS <- function(x, y, method, exact, conf.level){
 
-  x.var <- eval(parse(text=x), envir = parent.frame())
-  y.var <- eval(parse(text=y), envir = parent.frame())
+  x.var <- .loadServersideObject(x)
+  .checkClass(obj = x.var, obj_name = x, permitted_classes = c("numeric", "integer"))
+  y.var <- .loadServersideObject(y)
+  .checkClass(obj = y.var, obj_name = y, permitted_classes = c("numeric", "integer"))
   
   # get the number of pairwise complete cases
   n <- sum(stats::complete.cases(x.var, y.var))
@@ -26,9 +30,9 @@ corTestDS <- function(x, y, method, exact, conf.level){
   # runs a two-sided correlation test
   corTest <- stats::cor.test(x=x.var, y=y.var, method=method, exact=exact, conf.level=conf.level)
 
-  out <- list(n, corTest)
-  names(out) <- c("Number of pairwise complete cases", "Correlation test")
-  
+  out <- list(n, corTest, class = class(x.var))
+  names(out)[1:2] <- c("Number of pairwise complete cases", "Correlation test")
+
   # return the results
   return(out)
 
diff --git a/R/covDS.R b/R/covDS.R
index 9f645b62..fe25eabd 100644
--- a/R/covDS.R
+++ b/R/covDS.R
@@ -21,8 +21,10 @@
 #' of variables is not bigger than a percentage of the individual-level records (the allowed percentage is pre-specified
 #' by the 'nfilter.glm'). The second disclosure control checks that none of them is dichotomous with a level having fewer
 #' counts than the pre-specified 'nfilter.tab' threshold. If any of the input variables do not pass the disclosure
-#' controls then all the output values are replaced with NAs.
+#' controls then all the output values are replaced with NAs. The list also includes \code{class}, the class
+#' of the input object for client-side consistency checking.
 #' @author Amadou Gaye, Paul Burton, and Demetris Avraam for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #'
 covDS <- function(x=NULL, y=NULL, use=NULL){
@@ -36,14 +38,21 @@ covDS <- function(x=NULL, y=NULL, use=NULL){
   #nfilter.string <- as.numeric(thr$nfilter.string)
   #############################################################
   
-  x.val <- eval(parse(text=x), envir = parent.frame())
+  x.val <- .loadServersideObject(x)
+  .checkClass(obj = x.val, obj_name = x, permitted_classes = c("numeric", "integer", "matrix", "data.frame"))
+
   if (!is.null(y)){
-    y.val <- eval(parse(text=y), envir = parent.frame())
+    y.val <- .loadServersideObject(y)
+    .checkClass(obj = y.val, obj_name = y, permitted_classes = c("numeric", "integer", "matrix", "data.frame"))
   }
   else{
     y.val <- NULL
   }
-  
+
+  if (is.null(y.val) && any(class(x.val) %in% c("numeric", "integer"))) {
+    stop("If x is a numeric vector, y must also be a numeric vector.", call. = FALSE)
+  }
+
   # create a data frame for the variables
   if (is.null(y.val)){
     dataframe <- as.data.frame(x.val)
@@ -298,7 +307,7 @@ covDS <- function(x=NULL, y=NULL, use=NULL){
   
   }
 
-  return(list(sums.of.products=sums.of.products, sums=sums, complete.counts=complete.counts, na.counts=na.counts, errorMessage=errorMessage))
+  return(list(sums.of.products=sums.of.products, sums=sums, complete.counts=complete.counts, na.counts=na.counts, errorMessage=errorMessage, class=class(x.val)))
 
 }
 # AGGREGATE FUNCTION
diff --git a/R/kurtosisDS1.R b/R/kurtosisDS1.R
index 4f3f4e52..a51caa43 100644
--- a/R/kurtosisDS1.R
+++ b/R/kurtosisDS1.R
@@ -7,8 +7,9 @@
 #' @param method an integer between 1 and 3 selecting one of the algorithms for computing kurtosis
 #' detailed in the headers of the client-side \code{ds.kurtosis} function.
 #' @return a list including the kurtosis of the input numeric variable, the number of valid observations and
-#' the study-side validity message.
+#' \code{class}, the class of the input object for client-side consistency checking.
 #' @author Demetris Avraam, for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #'
 kurtosisDS1 <- function (x, method){
@@ -19,8 +20,9 @@ kurtosisDS1 <- function (x, method){
   nfilter.tab <- as.numeric(thr$nfilter.tab)
   #############################################################
   
-  x <- eval(parse(text=x), envir = parent.frame())
-  x <- x[stats::complete.cases(x)]
+  x.val <- .loadServersideObject(x)
+  .checkClass(obj = x.val, obj_name = x, permitted_classes = c("numeric", "integer"))
+  x <- x.val[stats::complete.cases(x.val)]
   
   if(length(x) < nfilter.tab){
     kurtosis.out <- NA
@@ -32,19 +34,16 @@ kurtosisDS1 <- function (x, method){
     
     if(method==1){
       kurtosis.out <- g2
-      studysideMessage <- "VALID ANALYSIS"
     }
     if(method==2){
       kurtosis.out <- ((length(x) + 1) * g2 + 6) * (length(x) - 1)/((length(x) - 2) * (length(x) - 3))
-      studysideMessage <- "VALID ANALYSIS"
     }
     if(method==3){
       kurtosis.out <- (g2 + 3) * (1 - 1/length(x))^2 - 3
-      studysideMessage <- "VALID ANALYSIS"
     }
   }
-  
-  out.obj <- list(Kurtosis=kurtosis.out, Nvalid=length(x), ValidityMessage=studysideMessage)
+
+  out.obj <- list(Kurtosis=kurtosis.out, Nvalid=length(x), class=class(x.val))
   return(out.obj)
   
 }
diff --git a/R/kurtosisDS2.R b/R/kurtosisDS2.R
index 1d4e3fec..e7790329 100644
--- a/R/kurtosisDS2.R
+++ b/R/kurtosisDS2.R
@@ -9,10 +9,10 @@
 #' @param global.mean a numeric, the combined mean of the input variable across all studies.
 #' @return a list including the sum of quartic differences between the values of x and the global mean of x across
 #' all studies, the sum of squared differences between the values of x and the global mean of x across all studies,
-#' the number of valid observations (i.e. the length of x after excluding missing values), and a validity message 
-#' indicating indicating a valid analysis if the number of valid observations are above the protection filter 
-#' nfilter.tab or invalid analysis otherwise.
+#' the number of valid observations (i.e. the length of x after excluding missing values), and \code{class},
+#' the class of the input object for client-side consistency checking.
 #' @author Demetris Avraam, for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #' 
 kurtosisDS2 <- function(x, global.mean){
@@ -23,20 +23,18 @@ kurtosisDS2 <- function(x, global.mean){
   nfilter.tab <- as.numeric(thr$nfilter.tab)
   #############################################################
   
-  x <- eval(parse(text=x), envir = parent.frame())
-  x <- x[stats::complete.cases(x)]
+  x.val <- .loadServersideObject(x)
+  .checkClass(obj = x.val, obj_name = x, permitted_classes = c("numeric", "integer"))
+  x <- x.val[stats::complete.cases(x.val)]
   
   if(length(x) < nfilter.tab){
-    sum_quartics.out <- NA
-    sum_squares.out <- NA
-    studysideMessage <- "FAILED: Nvalid less than nfilter.tab"
-  }else{
-    sum_quartics.out <- sum((x - global.mean)^4)
-    sum_squares.out <- sum((x - global.mean)^2)
-    studysideMessage <- "VALID ANALYSIS"
+    stop("FAILED: Nvalid less than nfilter.tab", call. = FALSE)
   }
-  
-  out.obj <- list(Sum.quartics=sum_quartics.out, Sum.squares=sum_squares.out, Nvalid=length(x), ValidityMessage=studysideMessage)
+
+  sum_quartics.out <- sum((x - global.mean)^4)
+  sum_squares.out <- sum((x - global.mean)^2)
+
+  out.obj <- list(Sum.quartics=sum_quartics.out, Sum.squares=sum_squares.out, Nvalid=length(x), class=class(x.val))
   return(out.obj)
   
 }
diff --git a/R/meanDS.R b/R/meanDS.R
index 59d1bc4e..576ec1b6 100644
--- a/R/meanDS.R
+++ b/R/meanDS.R
@@ -3,12 +3,15 @@
 #' @description Calculates the mean value.
 #' @details if the length of input vector is less than the set filter
 #' a missing value is returned.
-#' @param xvect a vector
-#' @return a numeric, the statistical mean
+#' @param x a character string, the name of a numeric or integer vector
+#' @return a list, with the estimated mean, the number of missing values, the number of
+#' valid values, the total number of values, and \code{class}, the class of the input
+#' object for client-side consistency checking
 #' @author Gaye A, Burton PR
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #'
-meanDS <- function(xvect){
+meanDS <- function(x){
 
   #############################################################
   # MODULE 1: CAPTURE THE nfilter SETTINGS
@@ -19,18 +22,18 @@ meanDS <- function(xvect){
   #nfilter.string <- as.numeric(thr$nfilter.string)
   #############################################################
 
+  xvect <- .loadServersideObject(x)
+  .checkClass(obj = xvect, obj_name = x, permitted_classes = c("numeric", "integer"))
+
   out.mean <- mean(xvect, na.rm=TRUE)
   out.numNa <- length(which(is.na(xvect)))
   out.totN <- length(xvect)
   out.validN <- out.totN-out.numNa
-  studysideMessage <- "VALID ANALYSIS"
-
   if((out.validN != 0) && (out.validN < nfilter.tab)){
-    out.mean <- NA
     stop("FAILED: Nvalid less than nfilter.tab", call. = FALSE)
   }
 
-  out.obj <- list(EstimatedMean=out.mean,Nmissing=out.numNa,Nvalid=out.validN,Ntotal=out.totN,ValidityMessage=studysideMessage)
+  out.obj <- list(EstimatedMean=out.mean,Nmissing=out.numNa,Nvalid=out.validN,Ntotal=out.totN,class=class(xvect))
   return(out.obj)
 
 }
diff --git a/R/meanSdGpDS.R b/R/meanSdGpDS.R
index 41fdb721..eecf0dde 100644
--- a/R/meanSdGpDS.R
+++ b/R/meanSdGpDS.R
@@ -3,17 +3,18 @@
 #' @description Server-side function called by ds.meanSdGp 
 #' @details Computes the mean and standard deviation across groups defined by one
 #' factor
-#' @param X a client-side supplied character string identifying the variable for which
+#' @param x a client-side supplied character string identifying the variable for which
 #' means/SDs are to be calculated
-#' @param INDEX a client-side supplied character string identifying the factor across
+#' @param index a client-side supplied character string identifying the factor across
 #' which means/SDs are to be calculated
 #' @author Burton PR
-#' 
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
+#'
 #' @return List with results from the group statistics
 #' @export
 #'
-meanSdGpDS <- function (X, INDEX){
-  
+meanSdGpDS <- function (x, index){
+
   #############################################################
   # MODULE 1: CAPTURE THE nfilter SETTINGS
   thr <- dsBase::listDisclosureSettingsDS()
@@ -23,9 +24,16 @@ meanSdGpDS <- function (X, INDEX){
   #nfilter.string <- as.numeric(thr$nfilter.string)
   #############################################################
 
+  X <- .loadServersideObject(x)
+  .checkClass(obj = X, obj_name = x, permitted_classes = c("numeric", "integer"))
+  INDEX <- .loadServersideObject(index)
+  .checkClass(obj = INDEX, obj_name = index, permitted_classes = c("factor", "character", "integer"))
+  x.class <- class(X)
+  index.class <- class(INDEX)
+
   FUN.mean <- function(x) {mean(x,na.rm=TRUE)}
   FUN.var <- function(x)  {stats::var(x,na.rm=TRUE)}
-  
+
   #Strip missings from both X and INDEX
   analysis.matrix<-cbind(X,INDEX)
   
@@ -114,8 +122,8 @@ meanSdGpDS <- function (X, INDEX){
   {
     table.valid<-TRUE
     cell.count.warning<-paste0("All tables valid") 
-    result<-list(table.valid,ansmat.mean,ansmat.sd,ansmat.count,Nvalid,Nmissing,Ntotal,cell.count.warning)
-    names(result)<-list("Table_valid","Mean_gp","StDev_gp", "N_gp","Nvalid","Nmissing","Ntotal","Message")
+    result<-list(table.valid,ansmat.mean,ansmat.sd,ansmat.count,Nvalid,Nmissing,Ntotal,cell.count.warning,x.class,index.class)
+    names(result)<-list("Table_valid","Mean_gp","StDev_gp", "N_gp","Nvalid","Nmissing","Ntotal","Message","class.x","class.index")
     return(result)
   }
   
@@ -123,8 +131,8 @@ meanSdGpDS <- function (X, INDEX){
   {
     table.valid<-FALSE
     cell.count.warning<-paste0("At least one group has between 1 and ", nfilter.tab-1, " observations. Please change groups") 
-    result<-list(table.valid,Nvalid,Nmissing,Ntotal,cell.count.warning)
-    names(result)<-list("Table_valid","Nvalid","Nmissing","Ntotal","Warning")
+    result<-list(table.valid,Nvalid,Nmissing,Ntotal,cell.count.warning,x.class,index.class)
+    names(result)<-list("Table_valid","Nvalid","Nmissing","Ntotal","Warning","class.x","class.index")
     return(result)
   }
   
diff --git a/R/quantileMeanDS.R b/R/quantileMeanDS.R
index 79fe3a96..b6a54e91 100644
--- a/R/quantileMeanDS.R
+++ b/R/quantileMeanDS.R
@@ -2,25 +2,31 @@
 #' @title Generates quantiles and mean information without maximum and minimum
 #' @description the probabilities 5%, 10%, 25%, 50%, 75%, 90%, 95% and the mean
 #' are used to compute the corresponding quantiles.
-#' @param xvect a numerical vector
-#' @return a numeric vector that represents the sample quantiles 
+#' @param x a character string, the name of a numeric or integer vector
+#' @return a list, with \code{quantiles}, a numeric vector that represents the sample
+#' quantiles, and \code{class}, the class of the input object for client-side consistency
+#' checking
 #' @export
 #' @author Burton, P.; Gaye, A.
-#' 
-quantileMeanDS <- function (xvect) {
-  
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
+#'
+quantileMeanDS <- function (x) {
+
+  xvect <- .loadServersideObject(x)
+  .checkClass(obj = xvect, obj_name = x, permitted_classes = c("numeric", "integer"))
+
   # check if the input vector is valid (i.e. meets DataSHIELD criteria)
   check <- isValidDS(xvect)
-  
+
   if(check){
-    # if the input vector is valid 
+    # if the input vector is valid
     qq <- stats::quantile(xvect,c(0.05,0.1,0.25,0.5,0.75,0.9,0.95), na.rm=TRUE)
     mm <- mean(xvect,na.rm=TRUE)
     quantile.obj <- c(qq, mm)
-    names(quantile.obj) <- c("5%","10%","25%","50%","75%","90%","95%","Mean")    
+    names(quantile.obj) <- c("5%","10%","25%","50%","75%","90%","95%","Mean")
   }else{
     quantile.obj <- NA
   }
-  
-  return(quantile.obj)
+
+  return(list(quantiles = quantile.obj, class = class(xvect)))
 }
diff --git a/R/skewnessDS1.R b/R/skewnessDS1.R
index 19f95dfc..64375375 100644
--- a/R/skewnessDS1.R
+++ b/R/skewnessDS1.R
@@ -7,8 +7,9 @@
 #' @param method an integer between 1 and 3 selecting one of the algorithms for computing skewness
 #' detailed in the headers of the client-side \code{ds.skewness} function.
 #' @return a list including the skewness of the input numeric variable, the number of valid observations and
-#' the study-side validity message.
+#' \code{class}, the class of the input object for client-side consistency checking.
 #' @author Demetris Avraam, for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #' 
 skewnessDS1 <- function(x, method){
@@ -19,8 +20,9 @@ skewnessDS1 <- function(x, method){
   nfilter.tab <- as.numeric(thr$nfilter.tab)
   #############################################################
   
-  x <- eval(parse(text=x), envir = parent.frame())
-  x <- x[stats::complete.cases(x)]
+  x.val <- .loadServersideObject(x)
+  .checkClass(obj = x.val, obj_name = x, permitted_classes = c("numeric", "integer"))
+  x <- x.val[stats::complete.cases(x.val)]
   
   if(length(x) < nfilter.tab){
     skewness.out <- NA
@@ -32,19 +34,16 @@ skewnessDS1 <- function(x, method){
     
     if(method==1){
       skewness.out <- g1
-      studysideMessage <- "VALID ANALYSIS"
     }
     if(method==2){
       skewness.out <- g1 * sqrt(length(x)*(length(x)-1))/(length(x)-2)
-      studysideMessage <- "VALID ANALYSIS"
     }
     if(method==3){
       skewness.out <- g1 * ((length(x)-1)/(length(x)))^(3/2)
-      studysideMessage <- "VALID ANALYSIS"
     }
   }
 
-  out.obj <- list(Skewness=skewness.out, Nvalid=length(x), ValidityMessage=studysideMessage)
+  out.obj <- list(Skewness=skewness.out, Nvalid=length(x), class=class(x.val))
   return(out.obj)
   
 }
diff --git a/R/skewnessDS2.R b/R/skewnessDS2.R
index 8d1cb484..aad3088d 100644
--- a/R/skewnessDS2.R
+++ b/R/skewnessDS2.R
@@ -9,10 +9,10 @@
 #' @param global.mean a numeric, the combined mean of the input variable across all studies.
 #' @return a list including the sum of cubed differences between the values of x and the global mean of x across
 #' all studies, the sum of squared differences between the values of x and the global mean of x across all studies,
-#' the number of valid observations (i.e. the length of x after excluding missing values), and a validity message 
-#' indicating indicating a valid analysis if the number of valid observations are above the protection filter 
-#' nfilter.tab or invalid analysis otherwise.
+#' the number of valid observations (i.e. the length of x after excluding missing values), and \code{class},
+#' the class of the input object for client-side consistency checking.
 #' @author Demetris Avraam, for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #' 
 skewnessDS2 <- function(x, global.mean){
@@ -23,21 +23,18 @@ skewnessDS2 <- function(x, global.mean){
   nfilter.tab <- as.numeric(thr$nfilter.tab)
   #############################################################
   
-  x <- eval(parse(text=x), envir = parent.frame())
-  x <- x[stats::complete.cases(x)]
+  x.val <- .loadServersideObject(x)
+  .checkClass(obj = x.val, obj_name = x, permitted_classes = c("numeric", "integer"))
+  x <- x.val[stats::complete.cases(x.val)]
   
   if(length(x) < nfilter.tab){
-    sum_cubes.out <- NA
-    sum_squares.out <- NA
-    studysideMessage <- "FAILED: Nvalid less than nfilter.tab"
-    stop(studysideMessage, call. = FALSE)
-  }else{
-    sum_cubes.out <- sum((x - global.mean)^3)
-    sum_squares.out <- sum((x - global.mean)^2)
-    studysideMessage <- "VALID ANALYSIS"
+    stop("FAILED: Nvalid less than nfilter.tab", call. = FALSE)
   }
-  
-  out.obj <- list(Sum.cubes=sum_cubes.out, Sum.squares=sum_squares.out, Nvalid=length(x), ValidityMessage=studysideMessage)
+
+  sum_cubes.out <- sum((x - global.mean)^3)
+  sum_squares.out <- sum((x - global.mean)^2)
+
+  out.obj <- list(Sum.cubes=sum_cubes.out, Sum.squares=sum_squares.out, Nvalid=length(x), class=class(x.val))
   return(out.obj)
   
 }
diff --git a/R/varDS.R b/R/varDS.R
index 390a9589..eb9fcf0e 100644
--- a/R/varDS.R
+++ b/R/varDS.R
@@ -3,15 +3,16 @@
 #' @description Calculates the variance.
 #' @details if the length of input vector is less than the set filter
 #' a missing value is returned.
-#' @param xvect a vector
+#' @param x a character string, the name of a numeric or integer vector
 #' @return a list, with the sum of the input variable, the sum of squares of the input variable,
 #' the number of missing values, the number of valid values, the number of total length of the
-#' variable, and a study message indicating whether the number of valid is less than the
-#' disclosure threshold
+#' variable, and \code{class}, the class of the input object for client-side
+#' consistency checking
 #' @author Amadou Gaye, Demetris Avraam, for DataSHIELD Development Team
+#' @author Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 #' @export
 #'
-varDS <- function(xvect){
+varDS <- function(x){
 
   #############################################################
   # MODULE 1: CAPTURE THE nfilter SETTINGS
@@ -22,21 +23,19 @@ varDS <- function(xvect){
   #nfilter.string <- as.numeric(thr$nfilter.string)
   #############################################################
 
+  xvect <- .loadServersideObject(x)
+  .checkClass(obj = xvect, obj_name = x, permitted_classes = c("numeric", "integer"))
+
   out.sum <- sum(xvect, na.rm=TRUE)
   out.sumSquares <- sum(xvect^2, na.rm=TRUE)
   out.numNa <- length(which(is.na(xvect)))
   out.totN <- length(xvect)
   out.validN <- out.totN-out.numNa
-  studysideMessage <- "VALID ANALYSIS"
-
   if((out.validN != 0) && (out.validN < nfilter.tab)){
-    out.sum <- NA
-    out.sumSquares <- NA
-    studysideMessage <- "FAILED: Nvalid less than nfilter.tab"
-    stop(studysideMessage, call. = FALSE)
+    stop("FAILED: Nvalid less than nfilter.tab", call. = FALSE)
   }
 
-  out.obj <- list(Sum=out.sum,SumOfSquares=out.sumSquares,Nmissing=out.numNa,Nvalid=out.validN,Ntotal=out.totN,ValidityMessage=studysideMessage)
+  out.obj <- list(Sum=out.sum,SumOfSquares=out.sumSquares,Nmissing=out.numNa,Nvalid=out.validN,Ntotal=out.totN,class=class(xvect))
   return(out.obj)
 
 }
diff --git a/man/corDS.Rd b/man/corDS.Rd
index 91e0a36d..6a734304 100644
--- a/man/corDS.Rd
+++ b/man/corDS.Rd
@@ -20,7 +20,8 @@ separately (columnwise) and the number of missing values casewise, and a vector
 sum of squares of each variable. The first disclosure control checks that the number of variables is
 not bigger than a percentage of the individual-level records (the allowed percentage is pre-specified
 by the 'nfilter.glm'). The second disclosure control checks that none of them is dichotomous with a
-level having fewer counts than the pre-specified 'nfilter.tab' threshold.
+level having fewer counts than the pre-specified 'nfilter.tab' threshold. The list also includes
+\code{class}, the class of the input object for client-side consistency checking.
 }
 \description{
 This function computes the sum of each vector of variable and the sum of the products
@@ -32,4 +33,6 @@ variables
 }
 \author{
 Paul Burton, and Demetris Avraam for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/corTestDS.Rd b/man/corTestDS.Rd
index 83c8ecc0..9c6d5cd2 100644
--- a/man/corTestDS.Rd
+++ b/man/corTestDS.Rd
@@ -22,7 +22,8 @@ only used for the Pearson product moment correlation coefficient if there are at
 4 complete pairs of observations.}
 }
 \value{
-the results of the correlation test.
+a list with the results of the correlation test and \code{class}, the class of the
+input object for client-side consistency checking.
 }
 \description{
 This function is similar to R function \code{cor.test}.
@@ -32,4 +33,6 @@ The function runs a two-sided correlation test
 }
 \author{
 Demetris Avraam, for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/covDS.Rd b/man/covDS.Rd
index 25b7d527..4b067a16 100644
--- a/man/covDS.Rd
+++ b/man/covDS.Rd
@@ -28,7 +28,8 @@ whether or not the input variables pass the disclosure controls. The first discl
 of variables is not bigger than a percentage of the individual-level records (the allowed percentage is pre-specified
 by the 'nfilter.glm'). The second disclosure control checks that none of them is dichotomous with a level having fewer
 counts than the pre-specified 'nfilter.tab' threshold. If any of the input variables do not pass the disclosure
-controls then all the output values are replaced with NAs.
+controls then all the output values are replaced with NAs. The list also includes \code{class}, the class
+of the input object for client-side consistency checking.
 }
 \description{
 This function computes the sum of each vector of variable and the sum of the products of each two
@@ -40,4 +41,6 @@ variables
 }
 \author{
 Amadou Gaye, Paul Burton, and Demetris Avraam for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/kurtosisDS1.Rd b/man/kurtosisDS1.Rd
index a6029a3d..239a3a43 100644
--- a/man/kurtosisDS1.Rd
+++ b/man/kurtosisDS1.Rd
@@ -14,7 +14,7 @@ detailed in the headers of the client-side \code{ds.kurtosis} function.}
 }
 \value{
 a list including the kurtosis of the input numeric variable, the number of valid observations and
-the study-side validity message.
+\code{class}, the class of the input object for client-side consistency checking.
 }
 \description{
 This function calculates the kurtosis of a numeric variable for each study separately.
@@ -25,4 +25,6 @@ The method is specified by the argument \code{method} in the client-side \code{d
 }
 \author{
 Demetris Avraam, for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/kurtosisDS2.Rd b/man/kurtosisDS2.Rd
index a61f16c4..19c07f95 100644
--- a/man/kurtosisDS2.Rd
+++ b/man/kurtosisDS2.Rd
@@ -14,9 +14,8 @@ kurtosisDS2(x, global.mean)
 \value{
 a list including the sum of quartic differences between the values of x and the global mean of x across
 all studies, the sum of squared differences between the values of x and the global mean of x across all studies,
-the number of valid observations (i.e. the length of x after excluding missing values), and a validity message 
-indicating indicating a valid analysis if the number of valid observations are above the protection filter 
-nfilter.tab or invalid analysis otherwise.
+the number of valid observations (i.e. the length of x after excluding missing values), and \code{class},
+the class of the input object for client-side consistency checking.
 }
 \description{
 This function calculates summary statistics that are returned to the client-side and 
@@ -29,4 +28,6 @@ of x across all studies and the number of valid observations of the input variab
 }
 \author{
 Demetris Avraam, for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/meanDS.Rd b/man/meanDS.Rd
index 6802ad58..c92c0d51 100644
--- a/man/meanDS.Rd
+++ b/man/meanDS.Rd
@@ -4,13 +4,15 @@
 \alias{meanDS}
 \title{Computes statistical mean of a vector}
 \usage{
-meanDS(xvect)
+meanDS(x)
 }
 \arguments{
-\item{xvect}{a vector}
+\item{x}{a character string, the name of a numeric or integer vector}
 }
 \value{
-a numeric, the statistical mean
+a list, with the estimated mean, the number of missing values, the number of
+valid values, the total number of values, and \code{class}, the class of the input
+object for client-side consistency checking
 }
 \description{
 Calculates the mean value.
@@ -21,4 +23,6 @@ a missing value is returned.
 }
 \author{
 Gaye A, Burton PR
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/meanSdGpDS.Rd b/man/meanSdGpDS.Rd
index 0b7cc1d5..031fdb14 100644
--- a/man/meanSdGpDS.Rd
+++ b/man/meanSdGpDS.Rd
@@ -4,13 +4,13 @@
 \alias{meanSdGpDS}
 \title{MeanSdGpDS}
 \usage{
-meanSdGpDS(X, INDEX)
+meanSdGpDS(x, index)
 }
 \arguments{
-\item{X}{a client-side supplied character string identifying the variable for which
+\item{x}{a client-side supplied character string identifying the variable for which
 means/SDs are to be calculated}
 
-\item{INDEX}{a client-side supplied character string identifying the factor across
+\item{index}{a client-side supplied character string identifying the factor across
 which means/SDs are to be calculated}
 }
 \value{
@@ -25,4 +25,6 @@ factor
 }
 \author{
 Burton PR
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/quantileMeanDS.Rd b/man/quantileMeanDS.Rd
index 5781685d..2126ed58 100644
--- a/man/quantileMeanDS.Rd
+++ b/man/quantileMeanDS.Rd
@@ -4,13 +4,15 @@
 \alias{quantileMeanDS}
 \title{Generates quantiles and mean information without maximum and minimum}
 \usage{
-quantileMeanDS(xvect)
+quantileMeanDS(x)
 }
 \arguments{
-\item{xvect}{a numerical vector}
+\item{x}{a character string, the name of a numeric or integer vector}
 }
 \value{
-a numeric vector that represents the sample quantiles
+a list, with \code{quantiles}, a numeric vector that represents the sample
+quantiles, and \code{class}, the class of the input object for client-side consistency
+checking
 }
 \description{
 the probabilities 5%, 10%, 25%, 50%, 75%, 90%, 95% and the mean
@@ -18,4 +20,6 @@ are used to compute the corresponding quantiles.
 }
 \author{
 Burton, P.; Gaye, A.
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/skewnessDS1.Rd b/man/skewnessDS1.Rd
index 76f48fc0..d042c9ee 100644
--- a/man/skewnessDS1.Rd
+++ b/man/skewnessDS1.Rd
@@ -14,7 +14,7 @@ detailed in the headers of the client-side \code{ds.skewness} function.}
 }
 \value{
 a list including the skewness of the input numeric variable, the number of valid observations and
-the study-side validity message.
+\code{class}, the class of the input object for client-side consistency checking.
 }
 \description{
 This function calculates the skewness of a numeric variable for each study separately.
@@ -25,4 +25,6 @@ The method is specified by the argument \code{method} in the client-side \code{d
 }
 \author{
 Demetris Avraam, for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/skewnessDS2.Rd b/man/skewnessDS2.Rd
index 4537f001..9b5e7b9b 100644
--- a/man/skewnessDS2.Rd
+++ b/man/skewnessDS2.Rd
@@ -14,9 +14,8 @@ skewnessDS2(x, global.mean)
 \value{
 a list including the sum of cubed differences between the values of x and the global mean of x across
 all studies, the sum of squared differences between the values of x and the global mean of x across all studies,
-the number of valid observations (i.e. the length of x after excluding missing values), and a validity message 
-indicating indicating a valid analysis if the number of valid observations are above the protection filter 
-nfilter.tab or invalid analysis otherwise.
+the number of valid observations (i.e. the length of x after excluding missing values), and \code{class},
+the class of the input object for client-side consistency checking.
 }
 \description{
 This function calculates summary statistics that are returned to the client-side and 
@@ -29,4 +28,6 @@ of x across all studies and the number of valid observations of the input variab
 }
 \author{
 Demetris Avraam, for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/man/varDS.Rd b/man/varDS.Rd
index 1c485c9b..022923be 100644
--- a/man/varDS.Rd
+++ b/man/varDS.Rd
@@ -4,16 +4,16 @@
 \alias{varDS}
 \title{Computes the variance of vector}
 \usage{
-varDS(xvect)
+varDS(x)
 }
 \arguments{
-\item{xvect}{a vector}
+\item{x}{a character string, the name of a numeric or integer vector}
 }
 \value{
 a list, with the sum of the input variable, the sum of squares of the input variable,
 the number of missing values, the number of valid values, the number of total length of the
-variable, and a study message indicating whether the number of valid is less than the
-disclosure threshold
+variable, and \code{class}, the class of the input object for client-side
+consistency checking
 }
 \description{
 Calculates the variance.
@@ -24,4 +24,6 @@ a missing value is returned.
 }
 \author{
 Amadou Gaye, Demetris Avraam, for DataSHIELD Development Team
+
+Tim Cadman, Genomics Coordination Centre, UMCG, Netherlands
 }
diff --git a/tests/testthat/test-disc-meanDS.R b/tests/testthat/test-disc-meanDS.R
index 22864733..41e3d9f9 100644
--- a/tests/testthat/test-disc-meanDS.R
+++ b/tests/testthat/test-disc-meanDS.R
@@ -25,7 +25,7 @@ set.standard.disclosure.settings()
 test_that("numeric meanDS, with below nfilter.tab values", {
     input <- c(NA, NA, 2.0, NA, 4.0)
     
-    expect_error(meanDS(input), "FAILED: Nvalid less than nfilter.tab", fixed = TRUE)
+    expect_error(meanDS("input"), "FAILED: Nvalid less than nfilter.tab", fixed = TRUE)
 })
 
 #
diff --git a/tests/testthat/test-disc-varDS.R b/tests/testthat/test-disc-varDS.R
index 3b60a771..28c8983d 100644
--- a/tests/testthat/test-disc-varDS.R
+++ b/tests/testthat/test-disc-varDS.R
@@ -25,7 +25,7 @@ set.standard.disclosure.settings()
 test_that("numeric varDS, with below nfilter.tab values", {
     input <- c(NA, NA, 2.0, NA, 4.0)
 
-    expect_error(varDS(input), "FAILED: Nvalid less than nfilter.tab", fixed = TRUE)
+    expect_error(varDS("input"), "FAILED: Nvalid less than nfilter.tab", fixed = TRUE)
 })
 
 #
diff --git a/tests/testthat/test-perf-meanDS.R b/tests/testthat/test-perf-meanDS.R
index 59266cb2..648ff3d3 100644
--- a/tests/testthat/test-perf-meanDS.R
+++ b/tests/testthat/test-perf-meanDS.R
@@ -36,7 +36,7 @@ test_that("numeric meanDS - performance", {
     .current.time <- .start.time
 
     while (difftime(.current.time, .start.time, units = "secs")[[1]] < .durationSec) {
-        meanDS(input)
+        meanDS("input")
 
         .count <- .count + 1
         .current.time <- Sys.time()
@@ -71,7 +71,7 @@ test_that("numeric meanDS, with NA - performance", {
     .current.time <- .start.time
 
     while (difftime(.current.time, .start.time, units = "secs")[[1]] < .durationSec) {
-        meanDS(input)
+        meanDS("input")
 
         .count <- .count + 1
         .current.time <- Sys.time()
diff --git a/tests/testthat/test-perf-varDS.R b/tests/testthat/test-perf-varDS.R
index 10fff94a..7abe84f6 100644
--- a/tests/testthat/test-perf-varDS.R
+++ b/tests/testthat/test-perf-varDS.R
@@ -36,7 +36,7 @@ test_that("numeric varDS - performance", {
     .current.time <- .start.time
 
     while (difftime(.current.time, .start.time, units = "secs")[[1]] < .durationSec) {
-        varDS(input)
+        varDS("input")
 
         .count <- .count + 1
         .current.time <- Sys.time()
@@ -71,7 +71,7 @@ test_that("numeric varDS, with NA - performance", {
     .current.time <- .start.time
 
     while (difftime(.current.time, .start.time, units = "secs")[[1]] < .durationSec) {
-        varDS(input)
+        varDS("input")
 
         .count <- .count + 1
         .current.time <- Sys.time()
diff --git a/tests/testthat/test-smk-corDS.R b/tests/testthat/test-smk-corDS.R
index bdc3607c..e86b6e2c 100644
--- a/tests/testthat/test-smk-corDS.R
+++ b/tests/testthat/test-smk-corDS.R
@@ -378,7 +378,7 @@ test_that("simple corDS, casewise, full", {
     res <- corDS("x", "y")
 
     expect_equal(class(res), "list")
-    expect_length(res, 5)
+    expect_length(res, 6)
 
     if (base::getRversion() < '4.0.0')
     {
@@ -465,7 +465,7 @@ test_that("simple corDS, casewise, neg. full", {
     res <- corDS("x", "y")
 
     expect_equal(class(res), "list")
-    expect_length(res, 5)
+    expect_length(res, 6)
 
     if (base::getRversion() < '4.0.0')
     {
@@ -552,7 +552,7 @@ test_that("simple corDS, casewise, some", {
     res <- corDS("x", "y")
 
     expect_equal(class(res), "list")
-    expect_length(res, 5)
+    expect_length(res, 6)
 
     if (base::getRversion() < '4.0.0')
     {
@@ -640,7 +640,7 @@ test_that("simple corDS, casewise, some", {
     res <- corDS("x", "y")
 
     expect_equal(class(res), "list")
-    expect_length(res, 5)
+    expect_length(res, 6)
 
     if (base::getRversion() < '4.0.0')
     {
@@ -718,6 +718,16 @@ test_that("simple corDS, casewise, some", {
     expect_equal(res$sums.of.squares[4], 58.0)
 })
 
+test_that("corDS throws error when object does not exist", {
+    expect_error(corDS("nonexistent_x", "nonexistent_y"), regexp = "does not exist")
+})
+
+test_that("corDS throws error when object is of invalid type", {
+    bad_input <- list(a = 1:3, b = 4:6)
+    y <- c(1.0, 2.0, 3.0)
+    expect_error(corDS("bad_input", "y"), regexp = "must be of type")
+})
+
 #
 # Done
 #
diff --git a/tests/testthat/test-smk-corTestDS.R b/tests/testthat/test-smk-corTestDS.R
index b500a085..cb54b62f 100644
--- a/tests/testthat/test-smk-corTestDS.R
+++ b/tests/testthat/test-smk-corTestDS.R
@@ -29,7 +29,7 @@ test_that("simple corTestDS, full, without na, pearson", {
     res <- corTestDS("x", "y", "pearson", NULL, 0.95)
 
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
 
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -82,7 +82,7 @@ test_that("simple corTestDS, neg. full, without na, pearson", {
     res <- corTestDS("x", "y", "pearson", NULL, 0.95)
 
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -135,7 +135,7 @@ test_that("simple corTestDS, some, pearson, without na, pearson", {
     res <- corTestDS("x", "y", "pearson", NULL, 0.95)
 
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -189,7 +189,7 @@ test_that("simple corTestDS, some, with na, pearson", {
     res <- corTestDS("x", "y", "pearson", NULL, 0.95)
 
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 3)
@@ -238,7 +238,7 @@ test_that("simple corTestDS, full, without na, kendall", {
     res <- corTestDS("x", "y", "kendall", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -284,7 +284,7 @@ test_that("simple corTestDS, neg. full, without na, kendall", {
     res <- corTestDS("x", "y", "kendall", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -330,7 +330,7 @@ test_that("simple corTestDS, some, kendall, without na, kendall", {
     res <- corTestDS("x", "y", "kendall", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -377,7 +377,7 @@ test_that("simple corTestDS, some, with na, kendall", {
     res <- corTestDS("x", "y", "kendall", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 3)
@@ -424,7 +424,7 @@ test_that("simple corTestDS, full, without na, spearman", {
     res <- corTestDS("x", "y", "spearman", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -470,7 +470,7 @@ test_that("simple corTestDS, neg. full, without na, spearman", {
     res <- corTestDS("x", "y", "spearman", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -516,7 +516,7 @@ test_that("simple corTestDS, some, spearman, without na, spearman", {
     res <- corTestDS("x", "y", "spearman", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 8)
@@ -563,7 +563,7 @@ test_that("simple corTestDS, some, with na, spearman", {
     res <- corTestDS("x", "y", "spearman", NULL, 0.95)
     
     expect_equal(class(res), "list")
-    expect_length(res, 2)
+    expect_length(res, 3)
     
     expect_length(res$`Number of pairwise complete cases`, 1)
     expect_equal(res$`Number of pairwise complete cases`, 3)
@@ -602,6 +602,16 @@ test_that("simple corTestDS, some, with na, spearman", {
     expect_equal(res$`Correlation test`$data.name[[1]], "x.var and y.var")
 })
 
+test_that("corTestDS throws error when object does not exist", {
+    expect_error(corTestDS("nonexistent_x", "nonexistent_y", "pearson", NULL, 0.95), regexp = "does not exist")
+})
+
+test_that("corTestDS throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    y <- c(1.0, 2.0, 3.0)
+    expect_error(corTestDS("bad_input", "y", "pearson", NULL, 0.95), regexp = "must be of type numeric or integer")
+})
+
 #
 # Done
 #
diff --git a/tests/testthat/test-smk-covDS.R b/tests/testthat/test-smk-covDS.R
index ce731938..24352ac7 100644
--- a/tests/testthat/test-smk-covDS.R
+++ b/tests/testthat/test-smk-covDS.R
@@ -27,7 +27,7 @@ test_that("numeric covDS, casewise.complete", {
 
     res <- covDS("input$v1", "input$v2", "casewise.complete")
 
-    expect_length(res, 5)
+    expect_length(res, 6)
     expect_equal(class(res), "list")
 
     res.sums.of.products.class <- class(res$sums.of.products)
@@ -130,7 +130,7 @@ test_that("numeric covDS, pairwise.complete", {
 
     res <- covDS("input$v1", "input$v2", "pairwise.complete")
 
-    expect_length(res, 5)
+    expect_length(res, 6)
     expect_equal(class(res), "list")
 
     res.sums.of.products.class <- class(res$sums.of.products)
@@ -232,6 +232,16 @@ test_that("numeric covDS, pairwise.complete", {
     expect_true(is.na(res$errorMessage))
 })
 
+test_that("covDS throws error when object does not exist", {
+    expect_error(covDS("nonexistent_x", "nonexistent_y", "pairwise.complete"), regexp = "does not exist")
+})
+
+test_that("covDS throws error when object is of invalid type", {
+    bad_input <- list(a = 1:3, b = 4:6)
+    y <- c(1.0, 2.0, 3.0)
+    expect_error(covDS("bad_input", "y", "pairwise.complete"), regexp = "must be of type")
+})
+
 #
 # Done
 #
diff --git a/tests/testthat/test-smk-kurtosisDS1.R b/tests/testthat/test-smk-kurtosisDS1.R
index fe939107..3a13786d 100644
--- a/tests/testthat/test-smk-kurtosisDS1.R
+++ b/tests/testthat/test-smk-kurtosisDS1.R
@@ -33,8 +33,7 @@ test_that("simple kurtosisDS1, method 1", {
     expect_equal(res$Kurtosis, -0.458210, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid,9)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("kurtosisDS1::smk::method 2")
@@ -49,8 +48,7 @@ test_that("simple kurtosisDS1, method 2", {
     expect_equal(res$Kurtosis, 0.270076, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid,9)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("kurtosisDS1::smk::method 3")
@@ -65,8 +63,16 @@ test_that("simple kurtosisDS1, method 3", {
     expect_equal(res$Kurtosis, -0.991672, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid,9)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
+})
+
+test_that("kurtosisDS1 throws error when object does not exist", {
+    expect_error(kurtosisDS1("nonexistent_object", 1), regexp = "does not exist")
+})
+
+test_that("kurtosisDS1 throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(kurtosisDS1("bad_input", 1), regexp = "must be of type numeric or integer")
 })
 
 #
diff --git a/tests/testthat/test-smk-kurtosisDS2.R b/tests/testthat/test-smk-kurtosisDS2.R
index 8f122a6e..69a735a6 100644
--- a/tests/testthat/test-smk-kurtosisDS2.R
+++ b/tests/testthat/test-smk-kurtosisDS2.R
@@ -36,8 +36,16 @@ test_that("simple kurtosisDS2", {
     expect_equal(res$Sum.squares, 3.25, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
+})
+
+test_that("kurtosisDS2 throws error when object does not exist", {
+    expect_error(kurtosisDS2("nonexistent_object", 2.5), regexp = "does not exist")
+})
+
+test_that("kurtosisDS2 throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(kurtosisDS2("bad_input", 2.5), regexp = "must be of type numeric or integer")
 })
 
 #
diff --git a/tests/testthat/test-smk-meanDS.R b/tests/testthat/test-smk-meanDS.R
index e6d81a73..8bb47c69 100644
--- a/tests/testthat/test-smk-meanDS.R
+++ b/tests/testthat/test-smk-meanDS.R
@@ -25,7 +25,7 @@ set.standard.disclosure.settings()
 test_that("numeric meanDS", {
     input <- c(0.0, 1.0, 2.0, 3.0, 4.0)
 
-    res <- meanDS(input)
+    res <- meanDS("input")
 
     expect_length(res, 5)
     expect_equal(class(res), "list")
@@ -37,15 +37,14 @@ test_that("numeric meanDS", {
     expect_equal(res$Nvalid, 5)
     expect_equal(class(res$Ntotal), "integer")
     expect_equal(res$Ntotal, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("meanDS::smk::numeric with NA")
 test_that("numeric meanDS, with NA", {
     input <- c(0.0, NA, 2.0, NA, 4.0)
 
-    res <- meanDS(input)
+    res <- meanDS("input")
 
     expect_length(res, 5)
     expect_equal(class(res), "list")
@@ -57,15 +56,14 @@ test_that("numeric meanDS, with NA", {
     expect_equal(res$Nvalid, 3)
     expect_equal(class(res$Ntotal), "integer")
     expect_equal(res$Ntotal, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("meanDS::smk::numeric with all NA")
 test_that("numeric meanDS, with all NA", {
-    input <- c(NA, NA, NA, NA, NA)
-    
-    res <- meanDS(input)
+    input <- rep(NA_real_, 5)
+
+    res <- meanDS("input")
     
     expect_length(res, 5)
     expect_equal(class(res), "list")
@@ -77,8 +75,16 @@ test_that("numeric meanDS, with all NA", {
     expect_equal(res$Nvalid, 0)
     expect_equal(class(res$Ntotal), "integer")
     expect_equal(res$Ntotal, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
+})
+
+test_that("meanDS throws error when object does not exist", {
+    expect_error(meanDS("nonexistent_object"), regexp = "does not exist")
+})
+
+test_that("meanDS throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(meanDS("bad_input"), regexp = "must be of type numeric or integer")
 })
 
 #
diff --git a/tests/testthat/test-smk-meanSdGpDS.R b/tests/testthat/test-smk-meanSdGpDS.R
new file mode 100644
index 00000000..de3d4ecb
--- /dev/null
+++ b/tests/testthat/test-smk-meanSdGpDS.R
@@ -0,0 +1,62 @@
+#-------------------------------------------------------------------------------
+# Copyright (c) 2019-2022 University of Newcastle upon Tyne. All rights reserved.
+# Copyright (c) 2022-2025 Arjuna Technologies, Newcastle upon Tyne. All rights reserved.
+#
+# This program and the accompanying materials
+# are made available under the terms of the GNU Public License v3.0.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program.  If not, see <http://www.gnu.org/licenses/>.
+#-------------------------------------------------------------------------------
+
+#
+# Set up
+#
+
+# context("meanSdGpDS::smk::setup")
+
+set.standard.disclosure.settings()
+
+#
+# Tests
+#
+
+# context("meanSdGpDS::smk::numeric by factor")
+test_that("simple meanSdGpDS, numeric by factor", {
+    x_var <- c(1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0)
+    index_var <- as.factor(c("A", "A", "A", "A", "A", "B", "B", "B", "B", "B"))
+
+    res <- meanSdGpDS("x_var", "index_var")
+
+    expect_equal(class(res), "list")
+    expect_true(res$Table_valid)
+    expect_equal(res$Nvalid, 10)
+    expect_equal(res$Nmissing, 0)
+    expect_equal(res$Ntotal, 10)
+    expect_equal(as.numeric(res$Mean_gp)[1], 3.0)
+    expect_equal(as.numeric(res$Mean_gp)[2], 8.0)
+})
+
+test_that("meanSdGpDS throws error when X does not exist", {
+    index_var <- as.factor(c("A", "A", "B", "B"))
+    expect_error(meanSdGpDS("nonexistent_x", "index_var"), regexp = "does not exist")
+})
+
+test_that("meanSdGpDS throws error when INDEX does not exist", {
+    x_var <- c(1.0, 2.0, 3.0, 4.0)
+    expect_error(meanSdGpDS("x_var", "nonexistent_index"), regexp = "does not exist")
+})
+
+test_that("meanSdGpDS throws error when X is not numeric or integer", {
+    bad_x <- c("a", "b", "c", "d")
+    index_var <- as.factor(c("A", "A", "B", "B"))
+    expect_error(meanSdGpDS("bad_x", "index_var"), regexp = "must be of type numeric or integer")
+})
+
+#
+# Done
+#
+
+# context("meanSdGpDS::smk::shutdown")
+
+# context("meanSdGpDS::smk::done")
\ No newline at end of file
diff --git a/tests/testthat/test-smk-quantileMeanDS.R b/tests/testthat/test-smk-quantileMeanDS.R
index 33eb0c6f..d305c8de 100644
--- a/tests/testthat/test-smk-quantileMeanDS.R
+++ b/tests/testthat/test-smk-quantileMeanDS.R
@@ -23,20 +23,25 @@
 test_that("numeric quantileMeanDS", {
     input <- c(0.0, 1.0, 2.0, 3.0, 4.0)
 
-    res <- quantileMeanDS(input)
+    res <- quantileMeanDS("input")
 
-    expect_length(res, 8)
-    expect_equal(class(res), "numeric")
-    expect_equal(res[[1]], 0.2)
-    expect_equal(res[[2]], 0.4)
-    expect_equal(res[[3]], 1.0)
-    expect_equal(res[[4]], 2.0)
-    expect_equal(res[[5]], 3.0)
-    expect_equal(res[[6]], 3.6)
-    expect_equal(res[[7]], 3.8)
-    expect_equal(res[[8]], 2.0)
+    expect_equal(class(res), "list")
+    expect_equal(res$class, "numeric")
 
-    res.names <- names(res)
+    qq <- res$quantiles
+
+    expect_length(qq, 8)
+    expect_equal(class(qq), "numeric")
+    expect_equal(qq[[1]], 0.2)
+    expect_equal(qq[[2]], 0.4)
+    expect_equal(qq[[3]], 1.0)
+    expect_equal(qq[[4]], 2.0)
+    expect_equal(qq[[5]], 3.0)
+    expect_equal(qq[[6]], 3.6)
+    expect_equal(qq[[7]], 3.8)
+    expect_equal(qq[[8]], 2.0)
+
+    res.names <- names(qq)
 
     expect_length(res.names, 8)
     expect_equal(class(res.names), "character")
@@ -54,20 +59,25 @@ test_that("numeric quantileMeanDS", {
 test_that("numeric quantileMeanDS, with NA", {
     input <- c(0.0, NA, 2.0, NA, 4.0)
 
-    res <- quantileMeanDS(input)
+    res <- quantileMeanDS("input")
+
+    expect_equal(class(res), "list")
+    expect_equal(res$class, "numeric")
+
+    qq <- res$quantiles
 
-    expect_length(res, 8)
-    expect_equal(class(res), "numeric")
-    expect_equal(res[[1]], 0.2)
-    expect_equal(res[[2]], 0.4)
-    expect_equal(res[[3]], 1.0)
-    expect_equal(res[[4]], 2.0)
-    expect_equal(res[[5]], 3.0)
-    expect_equal(res[[6]], 3.6)
-    expect_equal(res[[7]], 3.8)
-    expect_equal(res[[8]], 2.0)
+    expect_length(qq, 8)
+    expect_equal(class(qq), "numeric")
+    expect_equal(qq[[1]], 0.2)
+    expect_equal(qq[[2]], 0.4)
+    expect_equal(qq[[3]], 1.0)
+    expect_equal(qq[[4]], 2.0)
+    expect_equal(qq[[5]], 3.0)
+    expect_equal(qq[[6]], 3.6)
+    expect_equal(qq[[7]], 3.8)
+    expect_equal(qq[[8]], 2.0)
 
-    res.names <- names(res)
+    res.names <- names(qq)
 
     expect_length(res.names, 8)
     expect_equal(class(res.names), "character")
@@ -81,6 +91,15 @@ test_that("numeric quantileMeanDS, with NA", {
     expect_equal(res.names[[8]], "Mean")
 })
 
+test_that("quantileMeanDS throws error when object does not exist", {
+    expect_error(quantileMeanDS("nonexistent_object"), regexp = "does not exist")
+})
+
+test_that("quantileMeanDS throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(quantileMeanDS("bad_input"), regexp = "must be of type numeric or integer")
+})
+
 #
 # Done
 #
diff --git a/tests/testthat/test-smk-skewnessDS1.R b/tests/testthat/test-smk-skewnessDS1.R
index 562c3f65..f5d3357b 100644
--- a/tests/testthat/test-smk-skewnessDS1.R
+++ b/tests/testthat/test-smk-skewnessDS1.R
@@ -33,8 +33,7 @@ test_that("simple skewnessDS1, method 1", {
     expect_equal(res$Skewness, 0.443147, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid,9)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("skewnessDS1::smk::method 2")
@@ -49,8 +48,7 @@ test_that("simple skewnessDS1, method 2", {
     expect_equal(res$Skewness, 0.537175, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid,9)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("skewnessDS1::smk::method 3")
@@ -65,8 +63,16 @@ test_that("simple skewnessDS1, method 3", {
     expect_equal(res$Skewness, 0.3713805, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid,9)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
+})
+
+test_that("skewnessDS1 throws error when object does not exist", {
+    expect_error(skewnessDS1("nonexistent_object", 1), regexp = "does not exist")
+})
+
+test_that("skewnessDS1 throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(skewnessDS1("bad_input", 1), regexp = "must be of type numeric or integer")
 })
 
 #
diff --git a/tests/testthat/test-smk-skewnessDS2.R b/tests/testthat/test-smk-skewnessDS2.R
index 9e59061d..9e2ccea4 100644
--- a/tests/testthat/test-smk-skewnessDS2.R
+++ b/tests/testthat/test-smk-skewnessDS2.R
@@ -36,8 +36,16 @@ test_that("simple skewnessDS2", {
     expect_equal(res$Sum.squares, 3.25, tolerance = 1e-6)
     expect_equal(class(res$Nvalid), "integer")
     expect_equal(res$Nvalid, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
+})
+
+test_that("skewnessDS2 throws error when object does not exist", {
+    expect_error(skewnessDS2("nonexistent_object", 2.5), regexp = "does not exist")
+})
+
+test_that("skewnessDS2 throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(skewnessDS2("bad_input", 2.5), regexp = "must be of type numeric or integer")
 })
 
 #
diff --git a/tests/testthat/test-smk-varDS.R b/tests/testthat/test-smk-varDS.R
index 517b8d8f..51eac8e2 100644
--- a/tests/testthat/test-smk-varDS.R
+++ b/tests/testthat/test-smk-varDS.R
@@ -25,7 +25,7 @@ set.standard.disclosure.settings()
 test_that("numeric varDS", {
     input <- c(0.0, 1.0, 2.0, 3.0, 4.0)
 
-    res <- varDS(input)
+    res <- varDS("input")
 
     expect_length(res, 6)
     expect_equal(class(res), "list")
@@ -39,15 +39,14 @@ test_that("numeric varDS", {
     expect_equal(res$Nvalid, 5)
     expect_equal(class(res$Ntotal), "integer")
     expect_equal(res$Ntotal, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("varDS::smk::numeric with NA")
 test_that("numeric varDS, with NA", {
     input <- c(0.0, NA, 2.0, NA, 4.0)
 
-    res <- varDS(input)
+    res <- varDS("input")
 
     expect_length(res, 6)
     expect_equal(class(res), "list")
@@ -61,19 +60,18 @@ test_that("numeric varDS, with NA", {
     expect_equal(res$Nvalid, 3)
     expect_equal(class(res$Ntotal), "integer")
     expect_equal(res$Ntotal, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
 })
 
 # context("varDS::smk::numeric with all NA")
 test_that("numeric varDS, with all NA", {
-    input <- c(NA, NA, NA, NA, NA)
-    
-    res <- varDS(input)
-    
+    input <- rep(NA_real_, 5)
+
+    res <- varDS("input")
+
     expect_length(res, 6)
     expect_equal(class(res), "list")
-    expect_equal(class(res$Sum), "integer")
+    expect_equal(class(res$Sum), "numeric")
     expect_equal(res$Sum, 0)
     expect_equal(class(res$SumOfSquares), "numeric")
     expect_equal(res$SumOfSquares, 0)
@@ -83,8 +81,16 @@ test_that("numeric varDS, with all NA", {
     expect_equal(res$Nvalid, 0)
     expect_equal(class(res$Ntotal), "integer")
     expect_equal(res$Ntotal, 5)
-    expect_equal(class(res$ValidityMessage), "character")
-    expect_equal(res$ValidityMessage, "VALID ANALYSIS")
+    expect_equal(res$class, "numeric")
+})
+
+test_that("varDS throws error when object does not exist", {
+    expect_error(varDS("nonexistent_object"), regexp = "does not exist")
+})
+
+test_that("varDS throws error when object is not numeric or integer", {
+    bad_input <- c("a", "b", "c")
+    expect_error(varDS("bad_input"), regexp = "must be of type numeric or integer")
 })
 
 #