diff --git a/src/model_fitting/abundance_model.R b/src/model_fitting/abundance_model.R index 0f7c2dc13853b8cebd7d4dcc79defa729e82dbdb..8766ba31bf6460f09cd6951c8c2d93793af35f9e 100644 --- a/src/model_fitting/abundance_model.R +++ b/src/model_fitting/abundance_model.R @@ -388,8 +388,7 @@ results_res <- foreach(i = 1:nrow(all_model_terms), paste("SE", model_terms, sep = "_"), "theta", "SE.theta", "AIC", "R2" )} else { - result <- as.data.frame(matrix(NA, ncol = 3 * - length(model_terms) + 6, + result <- as.data.frame(matrix(NA, ncol = 3 * length(model_terms) + 6, nrow = 1)) names(result) <- c("model", paste("coeff", model_terms, sep = "_"), paste("P",model_terms,sep = "_"), @@ -446,17 +445,28 @@ results_res <- foreach(i = 1:nrow(all_model_terms), result[ , "R2"] <- summary(comparison_lm)$r.squared # if we are excluding years, this is the test of predicted data vs observed data # for this year (with which the model wasn't fitted) - if (!is.na(exclude_year)){ - predictors_excluded_year_pred <- predictors_excluded_year - predictors_excluded_year_pred$offset_term <- 0 - prediction_transect_excluded_year <- predict.glm(res, - newdata = predictors_excluded_year, - type = "response") - cross_lm = lm(log(predictors_excluded_year$ou_dens+ 1) ~ - log(prediction_transect_excluded_year + 1)) - - result[ , "R2_cross"] <- summary(cross_lm)$r.squared - } + if (!is.na(exclude_year)){ + predictors_excluded_year_pred <- predictors_excluded_year + predictors_excluded_year_pred$offset_term <- 0 + prediction_transect_excluded_year <- predict.glm(res, + newdata = predictors_excluded_year, + type = "response") + cross_lm_year = lm(log(predictors_excluded_year$ou_dens+ 1) ~ + log(prediction_transect_excluded_year + 1)) + + result[ , "R2_cross"] <- summary(cross_lm_year)$r.squared + } + if (!is.na(exclude_grid)){ + predictors_excluded_grid_pred <- predictors_excluded_grid + predictors_excluded_grid_pred$offset_term <- 0 + prediction_transect_excluded_grid <- predict.glm(res, + newdata = predictors_excluded_grid, + type = "response") + cross_lm_grid = lm(log(predictors_excluded_grid$ou_dens+ 1) ~ + log(prediction_transect_excluded_grid + 1)) + + result[ , "R2_cross"] <- summary(cross_lm_grid)$r.squared + } return(result) }