gsucarrat / gets

R Package for General-to-Specific (GETS) modelling and Indicator Saturation (ISAT) methods
8 stars 5 forks source link

Removing `isat` dependence on `Sys.call()` when determining to do Outlier Distortion #60

Closed moritzpschwarz closed 2 years ago

moritzpschwarz commented 2 years ago

This aims to deal with the multiple issues introduced by relying on the Sys.call() argument.

At the moment, this is only used explicitly in the part determining whether to carry out Outlier Distortion.

Motivating Error

Here is an example of the error:

data(Nile)
test_fun <- function(x = test1){
  isat(Nile, sis=TRUE, iis=TRUE, plot=TRUE, t.pval=0.005)
}
test_fun()

Solution Approach

I now don't rely on the Sys.call() argument anymore. Rather I'm creating a new list object called isat.args that is subsequently saved in getsis$aux$args.

This means that we can convert this:

## Outlier Proportion and Distortion Test 
  # Proportion Test was previously executed in the print.isat function
  if(!is.null(getsis$call$iis) & isTRUE(eval(getsis$call$iis))){
    if(!any(eval(getsis$call$sis), eval(getsis$call$tis), 
##change suggested by J-bat (implemented by G-man 12 June 2022):
            ifelse(is.logical(eval(getsis$call$uis)) & isTRUE(eval(getsis$call$uis)), TRUE, FALSE),
            !identical(userEstArg$name, "ols"))){

      getsis$outlier.proportion.test <- outliertest(getsis)
      getsis$outlier.distortion.test <- distorttest(getsis)
      }
    }

To this:

if(isat.args$iis & !any(isat.args$sis, isat.args$tis, isat.args$uis.logical)){
    getsis$outlier.proportion.test <- outliertest(getsis)
    getsis$outlier.distortion.test <- distorttest(getsis)
  }

Essential is for me that @jkurle checks if this works with his package. I think you added this fix here: !identical(userEstArg$name, "ols") - do we need this again? Happy to add in that case.

I'm also fixing an error that was raised by F-bear in 2020. The error was like this:

##uis as data.frame:
##in an email 5/10-2020, F-bear reported an error produced by
##the following code:
set.seed(123)
mx <- data.frame(matrix(runif(500),ncol=5))
colnames(mx) <- paste("x", seq(1:5), sep="")
mx <- as.data.frame(mx)
#mx <- as.zoo(mx) #conversion to zoo solves the issue partially
yx <- 2*mx[,1] + rnorm(100, 0, 1)
is2 <- isat(yx, iis=TRUE, sis=FALSE, uis=mx, t.pval=0.01)
is2

This used to fail as a data.frame was not treated the same as a matrix - a data.frame used to be broken up by column and each column was turned into a list. Now a data.frame is treated like a matrix.

jkurle commented 2 years ago

Yep, my code snippet is still required because your proportion and distortion tests only apply in the OLS case. The same code snippet seems to still work (without it, I get more failures of my unit tests) but even when I add it, I still get failures of other unit tests. So if it's fine with you, I will check it in more detail before confirming this review. It's probably just because we now have the additional getsis$aux$args but I'd like to make sure that it's nothing serious.

moritzpschwarz commented 2 years ago

Yep, my code snippet is still required

Excellent, thanks for the quick review. I've now added it again.

So if it's fine with you, I will check it in more detail before confirming this review. It's probably just because we now have the additional getsis$aux$args but I'd like to make sure that it's nothing serious.

Of course, take your time!!

jkurle commented 2 years ago

For me, the code does break some of my unit tests, other than just having the new output element getsis$aux$args.

When I call isat() from within my ivgets::ivisat() function, I get an error in line 56 when the list isat.args is created. The error message says "object t.pval not found". This is a bit of a mystery for me. Stepping into the function call before the error occurs, I see that t.pval exists in the environment but as a promise. The reason is that I am not specifying the argument t.pval in isat.default() explicitly but by passing the argument t.pval from ivgets::ivisat() to isat. This should not be a problem. The promise should simply be evaluated once it is used. For example, in the current version of isat.default(), I do not get any error messages even though it also uses the promise t.pval in later function calls.

I have also modified isat.default() with a simple line of code t.pval before line 56. Strangely enough, this does not raise an error. To be clear, the way I modified isat.default() is as follows:

t.pval
isat.args <- list(
    mc = mc, 
    ar = ar,
# and so on
    max.block.size = max.block.size, 
    t.pval = t.pval,
    wald.pval = wald.pval, 
# and so on
)

The first line t.pval runs without error (and simply prints t.pval to the console) but the next line raises the error "object t.pval not found". I'm puzzled.

jkurle commented 2 years ago

Okay, found the problem! I now have new branches for the ivgets and robust2sls packages that are compatible with this branch of gets. So once gets is updated to this new version, I can release or at least offer new, compatible versions accordingly.