\documentclass{article}

\usepackage{hyperref}

%\VignetteIndexEntry{Using the 'trouBBlme4SolveR' package}
%\VignettePackage{trouBBlme4SolveR}

\newcommand{\code}[1]{\texttt{#1}}
\newcommand{\pkg}[1]{{\normalfont\fontseries{b}\selectfont #1}}

\title{Using the \texttt{trouBBlme4SolveR} package}
\author{Iago Gin\'e-V\'azquez}

\begin{document}
\SweaveOpts{concordance=TRUE}

\maketitle

In this vignette we show an introduction to the package \pkg{trouBBlme4SolveR} and some examples of how to use the \code{dwmw} function (whose name was motivated because of \emph{Dealing With Model Warnings}).

In 2014, Ben Bolker wrote the publication \href{https://rpubs.com/bbolker/lme4trouble1}{https://rpubs.com/bbolker/lme4trouble1}, with some hints to solve convergence warnings produced by the functions \code{lmer} and \code{glmer}. Along the past years, he also have answered several related questions on the \pkg{lme4} repository in Github and in the SO forums. He also treated these issues in the \href{http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html}{GLMM FAQ}, mainly in the section \href{http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#troubleshooting}{Troubleshooting}.  This package was inspired by these documents and by the \pkg{lme4} documentation pages \texttt{troubleshooting} and \texttt{convergence}. This is the reason to make a homage to Ben Bolker in the package name, being a ``SolveR for (4) \pkg{lme4} troubles'', making the ``troub-lme4-SolveR'' a ``BB [Ben Bolker]-troub-lme4-SolveR'', i.e., \pkg{trouBBlme4SolveR}.

Let's start by the same example explained by Ben Bolker in his 2014's publication. Scaling and updating the optimizer to avoid model failed to converge is automatic by means of \code{dwmw}. Beyond that, while the final model in the publication is yet singular, the output model by \code{dwmw} is not.

<<<>>=
library(lme4)
data("fly_parameters", package = "trouBBlme4SolveR")
df <- fly_parameters
df$SUR.ID <- factor(df$SUR.ID)
df$replicate <- factor(df$replicate)
Rdet <- cbind(df$ValidDetections,df$FalseDetections)
Unit <- factor(1:length(df$ValidDetections))

m1 <- glmer(Rdet ~ tm:Area + tm:c.distance +
	    c.distance:Area + c.tm.depth:Area +
	    c.receiver.depth:Area + c.temp:Area +
	    c.wind:Area +
	    c.tm.depth + c.receiver.depth +
	    c.temp +c.wind + tm + c.distance + Area +
	    replicate +
	    (1|SUR.ID) + (1|Day) + (1|Unit) ,
    data = df, family = binomial(link="logit"))
summary(m1)
numcols <- grep("^c\\.",names(df))
dfs <- df
dfs[,numcols] <- scale(dfs[,numcols])
m1_sc <- update(m1,data=dfs)
ss <- getME(m1_sc,c("theta","fixef"))
m3 <- update(m1_sc,start=ss,control=glmerControl(optimizer="bobyqa",
						 optCtrl=list(maxfun=2e5)))
summary(m3)

library(trouBBlme4SolveR)
m1_new <- dwmw(m1, scale = TRUE, max_message_iter = 3)
summary(m1_new)
@


Next is an example in the \pkg{lme4} documentation, which is singular. Our function desingularizes it.

<<>>=
if(requireNamespace("nlme")){
	data(Orthodont,package="nlme")
	Orthodont$nsex <- as.numeric(Orthodont$Sex=="Male")
	Orthodont$nsexage <- with(Orthodont, nsex*age)
	fmo <- lmer(distance ~ age + (age|Subject) + (0 + nsex|Subject) +
		    (0 + nsexage|Subject), data = Orthodont)
	# without warnings
	fmo_new <- dwmw(fmo)
}
@

<<eval=FALSE>>=
summary(fmo)
@

<<echo=FALSE>>=
tryCatch(summary(fmo), error = function(e) "fmo object does not exist. Package 'nlme' should be installed first in your system.")
@

<<eval=FALSE>>=
summary(fmo_new)
@

<<echo=FALSE>>=
tryCatch(summary(fmo_new), error = function(e) "fmo_new object does not exist. Package 'nlme' should be installed first in your system.")
@



\subsection*{Other examples}

\begin{itemize}
	\item SO question \href{https://stackoverflow.com/questions/60028673/lme4-error-boundary-singular-fit-see-issingular}{lme4 error: boundary (singular) fit: see ?isSingular}
\end{itemize}

<<>>=
data("plants", package = "trouBBlme4SolveR")
fit <- lmer(Weight ~ 1 + (1|Rep:PLANT), data = plants)
summary(fit)
fit_new <- dwmw(fit)
summary(fit_new)
@

In this case, as the package does not analyze the random effect of each of the factors in an interaction among them (\code{Rep} and \code{PLANT}), it does not try to update the formula including them separately (\code{(1|Rep)} or \code{(1|PLANT)}), which is the final answer in the SO question, but it removes random effect specified and outputs a simple linear model.


\begin{itemize}
	\item \pkg{lme4} issue \href{https://github.com/lme4/lme4/issues/618}{convergence issues with continuous variables in model} at Github.
\end{itemize}

In this example, scaling the continuous predictor \href{https://github.com/lme4/lme4/issues/618#issuecomment-768589586}{makes the large-eigenvalue warning go away}.

<<>>=
data("issue618", package = "trouBBlme4SolveR")
fit <- glmer(outcome_dead ~ AGE + (1|ZIP), family = binomial, data = issue618)
summary(fit)
fit_new <- dwmw(fit, scale = TRUE)
summary(fit_new)
@

The same with the \href{https://github.com/lme4/lme4/issues/618#issuecomment-768672276}{larger dataset}:

<<>>=
data("issue618large", package = "trouBBlme4SolveR")
fit <- glmer(outcome_dead ~ AGE + (1|ZIP), family = binomial, data = issue618large)
summary(fit)
fit_new <- dwmw(fit, scale = TRUE)
summary(fit_new)
@


\begin{itemize}
	\item Cross Validated question \href{https://stats.stackexchange.com/questions/575666}{lme4: glmer() warning messages with count data mixed-effects model and how to proceed with model fit}
\end{itemize}

The convergence issue posted is solved by means of updating the model start parameters:

<<>>=
data("treatments", package = "trouBBlme4SolveR")
glmm.1 <- glmer(total_no ~ week * treatment * fzone + (1|plot), data = treatments, family = poisson)
summary(glmm.1)
glmm.11 <- dwmw(glmm.1, verbose = TRUE)
summary(glmm.11)
@




\begin{itemize}
	\item \href{https://rpubs.com/jimsavage/scale_issues}{A bag of tips and tricks for dealing with scale issues}
\end{itemize}

In this publication, the author suggests removing the convergence failing through dividing the variable \code{price} by \code{1000}. Another option is scaling (standardizing) all the continuous predictors.

<<>>=
if(requireNamespace("ggplot2")){
	data("diamonds", package = "ggplot2")

	# Grab the priciest diamonds
	diamonds_subset <- diamonds[(nrow(diamonds)-10000):nrow(diamonds),]
	# Fit the model
	fit_1 <- lmer(carat ~ depth + table + price + x + y + z + (1 + price | cut), data = diamonds_subset)
	# Let's try dividing price by 1000
	fit_2 <- lmer(carat ~ depth + table + I(price/1000) + x + y + z + (1 + I(price/1000) | cut), data = diamonds_subset)

	fit_new <- dwmw(fit_1, scale = TRUE, verbose = TRUE)
}
@


<<eval=FALSE>>=
summary(fit_1)
@

<<echo=FALSE>>=
tryCatch(summary(fit_1), error = function(e) "fit_1 object does not exist. Package 'ggplot2' should be installed first in your system.")
@


<<eval=FALSE>>=
summary(fit_2)
@

<<echo=FALSE>>=
tryCatch(summary(fit_2), error = function(e) "fit_2 object does not exist. Package 'ggplot2' should be installed first in your system.")
@


<<eval=FALSE>>=
summary(fit_new)
@

<<echo=FALSE>>=
tryCatch(summary(fit_new), error = function(e) "fit_new object does not exist. Package 'ggplot2' should be installed first in your system.")
@


\begin{itemize}
	\item SO question \href{https://stackoverflow.com/questions/45187681/how-to-use-update-for-random-part-in-lmer}{how to use update() for random part in lmer()?}
\end{itemize}

Function \code{fstruction} updates  the formula of singular models according to a similar proceeding to which is explained in that SO question.

\subsection{Session info}

<<<>>=
sessionInfo()
@


\end{document}