1 Motivational example

1.1 Birth weight study

Study on risks factor for low birth weights from Baystate medical center in Springfield (MA, USA) during 1986 [Hosmer & Lemeshow (1989), Applied Logistic Regression]

low: low birth weight indicator coded as: 1 (low birth weight, i.e. $\leq 2.5kg$) and 0 (normal birth weight)
age: mother’s age (in years) at the time of birth
lwt: mother’s weight (in pounds) at the time of the last menstrual period
race: mother’s race coded as: 1 (white), 2 (black) and 3 (other)
smoke: mother’s smoking status coded as 1 (smoker) and 0 (non-smoker)
ptl: number of previous premature labours
ht: hypertension history indicator coded as: 1(yes) and 0 (no)
ui: uterine irritability indicator coded as: 1 (yes) and 0 (no)
ftv: number of physician visits during the first trimester
bwt: infant’s birth weight (in grams)

2 `R` primer

2.1 `R` 101 for M2 PHDS

launch Rstudio
open a new .Rmd file from the default template
compile it with the "Knit" button
discuss with your neighbor to understand each part of the script

2.2 Importing data

import the dataset in birthweight.txt (you can use the "Import Dataset" button from Rstudio
briefly describe the data (use nice table outputs in Rmarkdown and ggplot2 graphics)

2.3 R basics

program a function that add the logarithm of two numbers
use browser() to debug your function when evaluating at negative arguments

3 Statistical modeling

3.1 Why do we do (bio)statistics ?

Statistics: summarizing information from experimental observations and quantifying the associated uncertainty

Always start with the research/scientific question !

3.2 How do we do (bio)statistics

ideally we would start from our scientific question, figure out what data we need to answer the question
usually the question is not very clear & the data have already been collected

Statistical Inference: we use a simple Generative Probabilistic model that could have generated the observations (Machine Learning sometimes reject this paradigm – cf. L. Breiman)

3.3 The likelihood

The likelihood is a fundamental building block of Biostatistics:

its maximization yields estimators that have good asymptocic (when $n$ the number of observations increases indefinitely) properties – no bias and smallest variance,
it also allows to perform significance testing of model parameters in many settings.

3.4 likelihood interpretation

The likelihood function quantifies how likely it is that a given (set of) observation(s) has been generated by our hypothesized Generative Probabilistic model.

The likelihood function is equal to the model joint probability distribution computed for the observations, and thus only a function of the model parameters.

3.5 Maximum Likelihood Estimator (MLE)

The idea of the MLE is to to optimize the likelihood function given the observations, by finding the model parameters that would give these observations the highest probability of being generated under the model.

seems like a reasonable and intuitive idea !

3.6 Bayesian statistics

Reverend Thomas Bayes proposed an alternative framework for statistical inference (actually before the “frequentist” method). It also relies on a probabilistic model through the likelihood function, but has different philosophical grounds than the frequentist.

To be continued…

4 Computational statistics

Computational statistics have become essential in modern statistics, with always bigger data, and always more sophisticated approaches.

4.1 Maximizing the likelihood

Maximizing the likelihood can easily be done analytically for simple linear models.

However: non-linear likelihoods are hard (sometimes impossible) to optimize analytically !

$\Rightarrow$ numerical optimization

4.2 Newton-Raphson optimizer

An algorithm to find values for which a function is zero.

Applied to the derivative of the likelihood function, this will identify the MLE
given that the log-likelihood is a concave function

4.3 Why always use the $log$ ?

Generally, we maximize the log-likelihood instead of the likelihood.

same thing ($\log$ is a monotonically increasing function)
products $\prod_{i=1}^n$ (e.g. likelihoods) become sums $\sum_{i=1}^n$
sums are easier to compute & easier to derivate
more precision for small positives numbers (e.g. small probabilities)

4.4 $log-exponential$ trick

This is not taught often but it can come very handy if you are writing your own statistical/optimisation program:

\[ \log \sum_{i=1}^n e^{x_i} = c + \log \sum_{i=1}^n e^{x_i-c} \]

Disclaimer: not useful today

4.5 Newton-Raphson optimizer intuition

start from an initial point $x_0$
tangential linear approximation of $f$ in $f(x_0)$
get the next point $x_1$ as the solution
…

Repeat until $x_{n+1}\approx x_n$

About the linear approximation at the $n + 1$ step: - goes through $f(x_n)$ - has slope $f'(x_n)$

So it has the following equation: $y= f'(x_n)(x-x_n) + f(x_n)$. Thus we find $x_{n+1}$ by setting $y=0$, which gives us $x_{n+1} = x_n -\frac{f(x_n)}{f'(x_n)}$

4.6 Newton-Raphson optimizer example

#install.packages("animation")
library("animation")
newton.method(FUN = function(x) (x - 2)^2 - 1, init = 9.5, 
              rg = c(-1, 10), tol = 0.001, interact = FALSE, 
              col.lp = c("orange", "red3", "dodgerblue1"), 
              lwd=1.5)

5 Practicals

5.1 What is the prevalence of low birth weights in our data ?

propose a generative probabilistic model
define the parameter of interest
Write down the log-likelihood and its the first derivative (by hand)
Show that the maximum value (that nullify the derivative) of this likelihood is for $\displaystyle\widehat{\pi}_{MLE} = \dfrac{1}{n}\sum_{i=1}^ny_i$

5.2 Brute force numerical optimization

Now let’s use R to program a Newton-Raphson algorithm to maximise this likelihood numerically

write an R function that computes either the likelihood or the log-likelihood for a probabiliyt p, with two additional arguments: i) obs the vector of observations (that is bw_data$low by default) and ii) log a logical (that is FALSE by default)
Write two R functions of p that compute the first and the second derivatives of the log-likelihood respectivelly (each with obs=bw_data$low as an additional argument)
plot those 4 functions (the likelihood, log-likelihood, )
write a Newton-Raphson function with 5 arguments (the first deriative of the function to maximize, its second derivative, the initial starting point, the tolerance, the maximum number of iterations)
use all three functions to compute the MLE of the low birthweight prevalence

5.3 Robustness and statistical inference

Try several initial values, different tolerance, and experiment with different stopping rules – i.e. convergence criterion – (you can combine them). Compare the results. Does it always work ?
Considering the null hypothesis $H_0: \pi=0.5$, compute the p-value for: i) the Wald test, ii)the score test, and iii) the Likelihood Ratio Test respectively, and compare to the glm and anova functions output.
Comment on the AIC and the 95% confidence interval

# logit and expit are bijective functions used to go from the real space ]-inf, +inf[ to the probability space [0,1] 
logit <- function(p){
  log(p/(1-p))
}
expit <- function(x){
  exp(x)/(1+exp(x))
}

# logisitc regression with just an intercept
my_logistic_reg <- glm(low~1, data=bw_data, family = "binomial")
summary(my_logistic_reg)
expit(my_logistic_reg$coefficients)

5.4 Does the mother’s smoking status impacts the probability of low birth weight

BONUS: repeat the exercise to assess wether the mother’s smoking status impacts the probability of low birth weight

Basic principles in Biostatistics: likelihood and statistical thinking

Boris Hejblum

September 10^th, 2019

1 Motivational example

1.1 Birth weight study

2 `R` primer

2.1 `R` 101 for M2 PHDS

2.2 Importing data

2.3 R basics

3 Statistical modeling

3.1 Why do we do (bio)statistics ?

3.2 How do we do (bio)statistics

3.3 The likelihood

3.4 likelihood interpretation

3.5 Maximum Likelihood Estimator (MLE)

3.6 Bayesian statistics

4 Computational statistics

4.1 Maximizing the likelihood

4.2 Newton-Raphson optimizer

4.3 Why always use the \(log\) ?

4.4 \(log-exponential\) trick

4.5 Newton-Raphson optimizer intuition

4.6 Newton-Raphson optimizer example

5 Practicals

5.1 What is the prevalence of low birth weights in our data ?

5.2 Brute force numerical optimization

5.3 Robustness and statistical inference

5.4 Does the mother’s smoking status impacts the probability of low birth weight

Basic principles in Biostatistics: likelihood and statistical thinking

Boris Hejblum

September 10th, 2019

1 Motivational example

1.1 Birth weight study

2 R primer

2.1 R 101 for M2 PHDS

2.2 Importing data

2.3 R basics

3 Statistical modeling

3.1 Why do we do (bio)statistics ?

3.2 How do we do (bio)statistics

3.3 The likelihood

3.4 likelihood interpretation

3.5 Maximum Likelihood Estimator (MLE)

3.6 Bayesian statistics

4 Computational statistics

4.1 Maximizing the likelihood

4.2 Newton-Raphson optimizer

4.3 Why always use the \(log\) ?

4.4 \(log-exponential\) trick

4.5 Newton-Raphson optimizer intuition

4.6 Newton-Raphson optimizer example

5 Practicals

5.1 What is the prevalence of low birth weights in our data ?

5.2 Brute force numerical optimization

5.3 Robustness and statistical inference

5.4 Does the mother’s smoking status impacts the probability of low birth weight

September 10^th, 2019

2 `R` primer

2.1 `R` 101 for M2 PHDS