Probability Distributions

R uses a consistent naming pattern for all distributions:
  d{dist}() — Probability density/mass function (PDF/PMF)
  p{dist}() — Cumulative distribution function (CDF)
  q{dist}() — Quantile function (inverse CDF)
  r{dist}() — Random number generation

Example for Normal distribution:
  dnorm(x, mean, sd)  — Density at x
  pnorm(q, mean, sd)  — P(X ≤ q)
  qnorm(p, mean, sd)  — Value where P(X ≤ value) = p
  rnorm(n, mean, sd)  — n random normal values

Available distributions:
  norm    — Normal
  binom   — Binomial
  pois    — Poisson
  t       — Student's t
  f       — F distribution
  chisq   — Chi-squared
  exp     — Exponential
  unif    — Uniform
  gamma   — Gamma
  beta    — Beta

library(ggplot2)

# Normal distribution parameters
mu <- 75000   # Mean salary
sigma <- 15000  # Standard deviation

# PDF — Density at a specific value
dnorm(75000, mean=mu, sd=sigma)  # Maximum (at mean)
dnorm(90000, mean=mu, sd=sigma)  # Lower (away from mean)

# CDF — Cumulative probability
pnorm(90000, mean=mu, sd=sigma)  # P(salary ≤ 90000) = 0.841
pnorm(60000, mean=mu, sd=sigma)  # P(salary ≤ 60000) = 0.159

# P(salary between 60K and 90K)
pnorm(90000, mu, sigma) - pnorm(60000, mu, sigma)  # 0.683 (68-95-99.7 rule!)

# Quantiles
qnorm(0.90, mean=mu, sd=sigma)  # 90th percentile = $94,224
qnorm(0.95, mean=mu, sd=sigma)  # 95th percentile
qnorm(0.975, mean=0, sd=1)      # 1.96 (used in 95% CI)

# Random generation
set.seed(42)
sample_salaries <- rnorm(1000, mean=mu, sd=sigma)
mean(sample_salaries)  # ≈ 75000
sd(sample_salaries)    # ≈ 15000

# Visualize normal distribution
x <- seq(30000, 120000, length.out=300)
df_norm <- data.frame(x=x, y=dnorm(x, mu, sigma))

ggplot(df_norm, aes(x, y)) +
  geom_line(color="#1565C0", size=1.5) +
  geom_area(data=subset(df_norm, x >= 60000 & x <= 90000),
             aes(y=y), fill="#90CAF9", alpha=0.4) +
  geom_vline(xintercept=mu, color="red", linetype="dashed") +
  scale_x_continuous(labels=function(x) paste0("$",x/1000,"K")) +
  labs(title="Salary Distribution (Normal)",
       subtitle="Shaded area = P(60K < X < 90K) ≈ 68%",
       x="Salary", y="Density") +
  theme_minimal()

# Binomial: n trials, probability p of success
# Example: 15 customers, 30% purchase probability
n <- 15; p <- 0.30

# P(exactly 5 buy) = P(X = 5)
dbinom(5, size=n, prob=p)   # 0.2061

# P(at most 5 buy) = P(X ≤ 5)
pbinom(5, size=n, prob=p)   # 0.7216

# P(more than 5 buy) = P(X > 5)
1 - pbinom(5, size=n, prob=p)  # 0.2784

# Expected value and variance
cat("E[X] =", n * p, "\n")         # 4.5
cat("Var[X] =", n * p * (1-p), "\n")  # 3.15
cat("SD[X] =", sqrt(n * p * (1-p)), "\n")  # 1.77

# Random binomial sample
rbinom(10, size=n, prob=p)

# Visualize
x_binom <- 0:15
df_binom <- data.frame(x=x_binom, prob=dbinom(x_binom, n, p))
ggplot(df_binom, aes(x, prob)) +
  geom_col(aes(fill=x==5), show.legend=FALSE) +
  scale_fill_manual(values=c("FALSE"="#90CAF9","TRUE"="#1565C0")) +
  labs(title="Binomial Distribution (n=15, p=0.3)",
       subtitle="Blue bar: P(X=5) = 0.206",
       x="Number of Purchases", y="Probability") +
  theme_minimal()

# Poisson: rate λ events per time period (count data)
# Example: λ = 4 customer arrivals per hour
lambda <- 4

# P(exactly 6 arrive) = P(X = 6)
dpois(6, lambda=lambda)   # 0.1042

# P(at most 4 arrive)
ppois(4, lambda=lambda)   # 0.6289

# P(more than 6 arrive)
1 - ppois(6, lambda=lambda)  # 0.1107

# Random Poisson
set.seed(42)
hourly_arrivals <- rpois(100, lambda=lambda)  # 100 hours
mean(hourly_arrivals)  # ≈ 4  (mean = lambda)
var(hourly_arrivals)   # ≈ 4  (Poisson: mean = variance!)

cat("Mean=Variance for Poisson:", round(mean(hourly_arrivals),2),
    "≈", round(var(hourly_arrivals),2), "\n")