---
title: "Logit model"
author: "Cédric NOEL - Jang SHILTZ"
email: "cedric.noel@univ-lorraine.fr"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
bibliography:
  - biblio.bib
vignette: >
  %\VignetteIndexEntry{Logit model}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```
# Introduction

In this section we consider a latent variable $y^*_{it}$ such that
\begin{equation}
y^*_{it} = f(a_{it}; \beta_k, \delta_k)+\epsilon_{it} =\beta_k A_{it}+\delta_k W_t+\epsilon_{it}
\end{equation}
where $\epsilon_{it}\sim \mathcal{N}\left(0;\ \sigma_k\right)$,  $A_{it}=(1,a_{it},a_{it}^2,\cdots,a_{it}^{n_\beta-1})^t $, $W_t=(w_{i1},\cdots,w_{in_\delta})^t $, $\beta_k =(\beta_{k1},\cdots, \beta_{kn_\beta})$ and $\delta_k=(\delta_{k1},\cdots,\delta_{kn_\delta}) $.\\
It is then classically assumed that the binary variable $y_{it}=1$ if $y^*_{it}>0$ and $y_{it}=0$ if $y^*_{it}\leq 0$.\\

If $\epsilon_{it}$ is assumed to follow a logistic distribution, we obtain the logit model.
Let $\rho_{ikt}=P(Y_{it}=1|W_i=w_i,C_i=k)$ be the probability of $y_{it}=1$ given membership in group $k$.
\begin{equation}
\rho_{ikt}=\dfrac{e^{\beta_k A_{it}+\delta_k W_{it}}}{1+e^{\beta_k A_{it}+\delta_k W_{it}}}
\end{equation}

# Parameters

We use artificial data contained in the library. These are 500 longitudinal data points of 0 and 1 values generated by the parameters below:

+ 3 groups ;
+ group probabilities are $\pi_1=0.22$, $\pi_2=0.44$ and $\pi_3=0.34$ ;
+ period is 10 time points ;
+ we use 3 polynomial shapes:
	- degree 3 with $\beta_1=(-2.266,-0.109,0.14,-0.011)$ ;
	- degree 4 with $\beta_2=(1.291,0.811,-0.554,0.077,-0.003)$ ;
	- degree 0 with $\beta_{3}=(-1.99243)$.


The data contains the time-dependent variable $Y_i$ in `data[,2:11]`, the time variable in `data[,12:21]`, a time-dependent covariate in `data[,24:33]`, and a time-invariant covariate that influences the group membership probability in `data[,48:49]`.

+ `data[,2:11]` is a matrix with 0 and 1 values.
+ `data[,12:21]` is a matrix with time 1 to 10.
+ `data[,24:33]` is a time-dependent covariate that influences the shape of the trajectory. It is a matrix with 0 and 1 values, representing the presence or absence of a characteristic for the individual.
+ `data[,48:49]` is a matrix with 0 and 1 values, representing a time-invariant covariate.


```{r}
library(trajeR)

data(data_LOGIT)

matplot(
  t(data_LOGIT[, 12:21]),t(data_LOGIT[, 2:11]),
  pch = 1, type='b', col="black", lty=1,
  xlab = "Times", ylab = "Values", main = "Plot of the individual's trajectories")
```

We use the each method to fit the model. For all method, we specify the number of group of our model, ng=3, the degree of the polynomial shape of the trajectories. Here we choice a line parallel to abscissa axis, a cubic polynomial and a quadric polynomial. So degree is vector $(0,3,4)$.

We specify `hessian=TRUE` to ask the calculus of the hessian matrix.

`Itermax` is set to 300 to assure a good approximation of the parameters.

For the Likelihood method we call `trajeR`  with option `Method ="L"`.

```{r, message = FALSE}
solL <- trajeR(
  Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21],
  ng = 3, degre = c(0,3,4),
  Model = "LOGIT", Method = "L",
  hessian = TRUE
  )

solL
```

For the EM method, we use the same parameters except `Method = EM` or `EMIRLS`. We specifiy `itermax = 300` to ensure convergence.

```{r, message = FALSE}
#EM
solEM <- trajeR(Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21],
               ng = 3, degre = c(0,3,4), 
               Model = "LOGIT", Method = "EM", hessian = TRUE, itermax = 300) 
#EMIRLS
solEMIRLS <- trajeR(Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21],
                    ng = 3, degre=c(0,3,4),
                    Model = "LOGIT", Method = "EMIRLS", hessian = TRUE,
                    itermax = 300)

solEM
solEMIRLS
```
We can add risk covariate that influence the belonging probability. By default, the effect of risk variable are compared to the first group that is reference group but we can change this setting by the option `refgr`.



```{r}
solLRisk <- trajeR(Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21], Risk = data_LOGIT[,22:23],
                  ng = 3, degre = c(0,3,4),
                  Model = "LOGIT", Method = "L", hessian = TRUE,
                  itermax = 300)

solLRisk
```


# Plot

We plot only the trajectory on a graph. We have just to use the `trajeR` plot function. It use the class of the object to draw the polynomial shape. The values in abscissa are the first row of the time variable.

```{r}
plotrajeR(solLRisk)
```

We can add longitudinal data to this graph.
If we want them on the plot we have to specify Y and A in the function `plot()`.
For more visibility, we have enlarge the 0 and 1 value on the graph and plot the data with a little and random shift, control by the index `dec` in the function plot.

By default colors are gray scale, but we can specify colors we want.

```{r}
# colour's defintion
trans <- "70"
col1 <- "#034569"
col1.1 <- paste0("#64AAD0", trans)
col2 <- "#750062"
col2.1 <- paste0("#D962C7", trans)
col3 <- "#A68900"
col3.1 <- paste0("#FFE773", trans)
cols1 <- c(col1.1, col2.1, col3.1)
cols2 <- c(col1, col2, col3)
vcol <- c(cols1, cols2)

plotrajeR(solLRisk, Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21], dec = 5, col = vcol)
```

A time-dependent covariate can directly influence the shape of the trajectories. In the data, we have a time-dependent covariate composed of 0 and 1 values. We can account for its effects by using the `TCOV` option in the `trajeR` command. Any of the three methods above can be used.


```{r}
solLTCOV <- trajeR(
  Y = data_LOGIT[,2:11], A =data_LOGIT[,12:21], TCOV = data_LOGIT[,24:33],
  ng = 3, degre = c(0,3,4),
  Model = "LOGIT", Method = "L", hessian = TRUE,
  itermax = 300
  )

solEMTCOV <- trajeR(
  Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21], TCOV = data_LOGIT[,24:33],
  ng = 3, degre = c(0,3,4),
  Model = "LOGIT", Method = "EM", hessian = TRUE,
  itermax = 300
  )

solEMIRLSTCOV <- trajeR(
  Y = data_LOGIT[,2:11], A = data_LOGIT[,12:21], TCOV = data_LOGIT[,24:33],
  ng = 3, degre = c(0,3,4),
  Model = "LOGIT", Method = "EMIRLS", hessian = FALSE,
  itermax = 300
  )

solLTCOV
solEMTCOV
solEMIRLSTCOV

plotrajeR(solLTCOV, col = vcol)
```

If we want show the impact of a particular value of the time covariate in the trajectory, we can add this to the plot by `plotcov` option.
The fill line is the trajectory with the time covariate matrix and the dashed one show the impact on this trajectory of a particular value.


```{r}
plotrajeR(solLTCOV, col = vcol, plotcov = c(0,0,0,0,0,1,1,1,1,1))
```