Introduction to Statistical Inference

  • Statistical terms and key concepts.
    • population
    • sample
    • probability distributions
    • population parameters
    • point estimates
    • standard error
    • assumptions
  • Goal:
    • Understand how they are all related
    • Know how to work with them in R

Overview

    1. What is a population?
    2. What is a sample?
    3. How does the sample relate to the population?

1. What is a population?

In statistics, a population is the entire group of individuals, objects, measurements or outcomes that you want to learn about.

If you had the data for the entire population, you would be able ask and answer very direct questions.

  • Imagine every citizen of US had a chip implanted in their body that could instantaneosly measure their current blood pressure
    • How many people in united states have a systolic BP > 120?
    • What is the median blood pressure?
  • Complete observation of the entire population is usually not possible
    • It’s too expensive to collect everything
    • Some questions would be hard to answer even if we did have all the observations from the full population
  • Instead, we try to represent the distribution of a population mathematically
    • this is a simplification of the data, where we try to capture the essence of a collection of numbers with a much smaller set of numbers (parameters).

In the example below, we read in a dataset of 1000 numbers, draw a histogram of their values, then overlay a probability density function. The probability function is completely defined by just 2 parameters and closely matches the distribution of the dataset.

numdat = read.csv("numdata1k.csv")
dim(numdat)
## [1] 1000    1
summary(numdat)
##        x         
##  Min.   : 59.21  
##  1st Qu.: 91.85  
##  Median :100.44  
##  Mean   :100.62  
##  3rd Qu.:109.18  
##  Max.   :152.43
hist(numdat[,1], probability=T, breaks=100)
curve(dnorm(x, mean = 100, sd = 13), 
      col = "red", lwd = 2, add = TRUE)

If we can figure out a way to represent the data as a mathematical function, then we can ask questions using that mathematical function instead of having to deal with the full set of numbers.

In the plot above, the red line comes from a normal(Gaussian) probability distribution with mean=100 and standard dev=13. It’s a pretty good fit for the numbers we just read in.

A normal(Gaussian) distribution can be fully described by two parameters:

  • \(f(x) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left( -\frac{(x - \mu)^2}{2\sigma^2} \right)\)
  • \(\mu\) : the mean (center)
  • \(\sigma\) : standard deviation (spread)
  • The area under the curve always sums to 1.

If you know the 2 parameters \(\mu\) and \(\sigma\), then you know everything there is know about the distribution.

That red line captures the ‘essence’ of 1000 numbers with just 2 parameters.

.

Calculus Quiz 1

    • \(f(x) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left( -\frac{(x - \mu)^2}{2\sigma^2} \right)\)
    • \(\mu\) : the mean (center)
    • \(\sigma\) : standard deviation (spread)
  • Q: If your population data follows a normal distribution with mean \(\mu\) = 100 and standard deviation \(\sigma\)=13, what’s the probability of drawing a random observation greater than 120?
  • Show answer A: \(\int_{120}^{\infty} \frac{1}{\sqrt{2\pi 13^2}} \exp\left( -\frac{(x - 100)^2}{2*13^2} \right) dx\)

  • calculate this probability in R

     1 - pnorm (120, mean=100, sd=13)
    ## [1] 0.0619679

  • We can also use R to visualize the problem.

    # Parameters
mu <- 100
sigma <- 13

# Create a sequence of x values
x_vals <- seq(mu - 4*sigma, mu + 4*sigma, length = 500)

# Normal density values
y_vals <- dnorm(x_vals, mean = mu, sd = sigma)

# Base plot
plot(x_vals, y_vals, type = "l", lwd = 2, col = "black",
     ylab = "Density", xlab = "x", main = "P(X > 120) under Normal Curve")

# Shade area where x > 120
x_shade <- seq(120, max(x_vals), length = 100)
y_shade <- dnorm(x_shade, mean = mu, sd = sigma)
polygon(c(120, x_shade, max(x_shade)),
        c(0, y_shade, 0),
        col = "lightblue", border = NA)

# Add vertical line at x = 120
abline(v = 120, col = "red", lty = 2)

# Add legend with computed probability
p_gt_120 <- 1 - pnorm(120, mean = mu, sd = sigma)
legend("topright", legend = paste0("P(X > 120) = ", round(p_gt_120, 4)),
       bty = "n")

.

Why don’t we work with population?

    • In reality
      • We don’t have the full details of the population
      • The observations don’t really follow a known mathematical distribution
      • If they do follow a distribution, we don’t know the exact parameters
    • In practice
      • we use a subset of the population
      • We make assumptions about how we think the population is distributed
      • we use the subset to figure out details(parameters) of the population distribution
    • Can we draw meaningful conclusions about the population parameters?
      • Was the study designed properly?
      • is the sample big enough?
      • are our assumptions about the population correct?
      • are we using the right testing procedure?

2. What is a sample?

  • It’s a subset of population that we use to draw conclusions about the population.

  • In studying lung cancer in the Bronx, it’s a lot easier to collect info on 30 patients than every single patient that’s ever been treated here at Montefiore and the surrounding hospitals.

Point Estimates

  • how do we make inferences about the population parameters?
    • Points estimates computed from your sample.
  • What are the point estimates?
    • It all depends on the assumptions made of the population distribution and what you’re trying to determine
  • What are some examples of point estimates?
    • sample mean: \(\bar{x} = \frac{1}{n} \sum_{i=1}^{n} x_i\)
    • sample variance: \(s^2 = \frac{1}{n - 1} \sum_{i=1}^{n} (x_i - \bar{x})^2\)
    • sample proportion: \(\hat{p} = \frac{x}{n}\)

.

3. How does the sample relate to the population?

    • We use the sample point estimates to figure out the population parameters.

      • if the sample mean is 75, we would hope the population mean is also close to 75.

    • When we take a sample from a population and calculate a point estimate (like a sample mean or proportion), the result will vary from sample to sample because of the natural randomness in the sampling process.

      • if you repeated this procedure many times, you wouldn’t always get the same point estimate.

    • the standard error (SE) is a measurement of the spread in the sample point estimate.

    • Remember this distinction:

      • standard deviation (SD) is a general function that measures the spread of a set of numbers. It is usually used to describe a population parameter.
      • standard error (SE) is only relevant when you’re talking about sample point estimates and their relation to a population parameter. It’s the amount of spread you would expect to see from the point estimate if you repeated the sampling procedure many times. The larger the sample size, the less spread you’ll see and the closer to the true population parameter you’ll be (for unbiased estimators).

Statistical Inference

    • If you know the population follows a certain distribution
      • Then you can figure out the distribution of a sample point estimate.
    • If the population observations are independent and identically distributed and follow a normal distribution with parameters \(\mu\) and \(\sigma\),
      • Let \(x_1, x_2, \ldots, x_n \overset{\text{i.i.d.}}{\sim} \mathcal{N}(\mu, \sigma^2)\)
      • Then the sample mean also follows a normal distribution.
        • \(\bar{x} \sim \mathcal{N}\left(\mu, \frac{\sigma^2}{n} \right)\)
    • If you know the distribution of the point estimator
      • Then you can ask mathematical questions about the population parameters.
    • What are the usual questions?
      • What’s the probability that an interval within a certain range of the sample estimate includes the true population parameter?
        • confidence intervals
      • Assuming there’s no effect in the population, how likely is it to observe a sample point estimate that is as extreme(or more extreme) just by chance?
        • hypothesis tests and p-values
      • How big does my sample size need to be in order to be able to accurately say things about the population?
        • power calculations
    • logic statement (\(A \Rightarrow B\))
      • if A is true, then B is true
        • what if A is False?
    • If the statistical assumptions of our tests are true then we can draw conclusions about the population parameters based on the sample point estimator.
      • What if the assumptions are not met?
    • If \(x_1, x_2, \ldots, x_n \overset{\text{i.i.d.}}{\sim} \mathcal{N}(\mu, \sigma^2))\) then \(\bar{x} \sim \mathcal{N}\left(\mu, \frac{\sigma^2}{n} \right)\) .
      • What if the population isn’t normally distributed?
    • Our statistical conclusions are only valid if the assumptions of the test are met

.

How do you check that your assumptions on the population are correct?

Without having the complete observations of the entire population at hand it will be very hard to prove the underlying assumptions you make about the population parameters are true. The closest we can do is to show they are reasonable.

There are standard procedures for testing the assumptions of every common statistical test (e.g., Shapiro Wilk, Bartlett’s),

But it’s always a good idea to start with a visualization that shows how the distribution of your sample fairs against the distribution you’re assuming it to take.

QQ plot.

Stands for the quantile-quantile plot. Here we show a scatterplot of the quantiles of our sample with the expected quantiles of the distribution we’re expecting it to follow. If the plots coincide with each other, then we can proceed with our test. If not, we should figure out a different test that doesn’t violate its assumptions.

# Step 1: Simulate some data
x <- rnorm(100, mean = 5, sd = 2)  # Normally distributed sample

# Step 2: Sort the sample data
x_sorted <- sort(x)

# Step 3: Compute theoretical quantiles from the standard normal
n <- length(x_sorted)
# Use (i - 0.5)/n quantile levels
p <- (1:n - 0.5) / n
theoretical_q <- qnorm(p, mean = mean(x), sd = sd(x))  # match sample mean/sd

# Step 4: Plot
plot(theoretical_q, x_sorted,
     main = "QQ Plot from Scratch",
     xlab = "Theoretical Quantiles",
     ylab = "Sample Quantiles",
     pch = 19, col = "blue")

# Step 5: Add reference line (through first and third quartiles)
q_theo <- quantile(theoretical_q, c(0.25, 0.75))
q_sample <- quantile(x_sorted, c(0.25, 0.75))
slope <- diff(q_sample) / diff(q_theo)
intercept <- q_sample[1] - slope * q_theo[1]

abline(intercept, slope, col = "red", lwd = 2)

But in practice, you would just use the built-in functions. To check if your data follows a normal distribution, you can use the qqnorm, and shapiro.test functions.

#visual check
qqnorm(x)

#procedural check
#if the p-value is small, then your assumptions are violated
shapiro.test(x)
## 
##  Shapiro-Wilk normality test
## 
## data:  x
## W = 0.98778, p-value = 0.4921

How do you work with probability distributions in R?

In the past we had to use lookup tables to find our probabilities.

Now we use computers.

In R, the functions for working with probability distributions all have a similar format.

  • Normal (Gaussian) Distribution
    • dnorm(x, …): compute the probability density at x
    • pnorm(q, …): compute the cumulative probability density at q.
      The area under the curve to the left of point q.
    • qnorm(p, …): compute the inverse cumulative probability. what is the value of x that would give you an area under the curve of p?
    • rnorm(n, …): generate random numbers
  • format for calling probability functions
    • start with a prefix letter that tells what you want: d, p, q, r
    • end with a suffix for the name of the distribution

distributions

Distribution Name in R Notes
Normal norm Continuous
Binomial binom Discrete
Poisson pois Discrete
Exponential exp Continuous
Chi-squared chisq Continuous
Student’s t t Continuous
F-distribution f Continuous
Uniform unif Continuous
Geometric geom Discrete
Negative Binomial nbinom Discrete
Hypergeometric hyper Discrete
Gamma gamma Continuous
Beta beta Continuous
Weibull weibull Continuous
Logistic logis Continuous

examples

What’s the area under the curve to the the left of 0, in a standard normal distribution?

pnorm(0)
## [1] 0.5

What value of x would give you an area under the curve of 0.5?

qnorm(0.5)
## [1] 0

Plot the values from -3 to 3 of the standard normal distribution

#first get the range of x values you to to plot
xvals = seq(-3, 3, 0.01)

#then get the probabilites of each of those values
yvals = dnorm(xvals)

#now plot them in a scatter plot
plot(xvals, yvals)

Find the area to the right of 1.96 in a standard normal

#this is the area to the left
pnorm(1.96)
## [1] 0.9750021
#to get the area to the right, you can subtract it from 1
1 - (pnorm(1.96))
## [1] 0.0249979

What’s the area between -1.96 and 1.96 in a standard normal?

pnorm(1.96) - pnorm(-1.96)
## [1] 0.9500042

What’s the area under the curve of all the values higher than 1.96, or lower than -1.96?

#area to the left of the lower point
a1 = pnorm(-1.96)
a1
## [1] 0.0249979
#area to the right of the higher point
a2 = 1 - pnorm(1.96)
a2
## [1] 0.0249979
a1+a2
## [1] 0.04999579

generate 30 random numbers from a standard normal distribution

x = rnorm(30)
x
##  [1] -1.14918688  0.44046906 -0.20203362  0.45611542 -1.92099665 -2.70246439
##  [7] -0.37669252  0.38939779 -0.07782196 -1.00353317 -0.25170764  0.37935000
## [13] -1.08887781  1.07647682 -0.31622814  1.40120373  0.33669312  1.36952001
## [19]  0.57901692  0.06778218 -0.70083398 -0.27488586 -0.90526724 -0.10986025
## [25]  0.44286151 -0.79185606  0.43655710  0.30467749  1.26903804  0.29551744

.

What are the next the lectures about?

  • Finding the right statistical test for your data
  • checking the assumptions
  • how to run the test.
  • interpreting the results