effect size, P value, and Bayes odds

P value, effect size, and Bayes factor, and after one year

• 2021
• Cole (2020): Surprise! (P value vs. S value)
• KDnuggets by Bills (2021): Null Hypothesis Significance Testing is Still Useful
• 2019:Another wave beyond the P-value:
• Greenland (2019): Valid p-values behave exactly as they should: Some misleading criticisms of p-values and their resolution with s-values.
• Benjamin (2019): Three recommendations for improving the use of p-values
• Blume (2019): An introduction to second-generation p-values
• Tarran (2019): Is this the end of "statistical significance"?
• American Statistician (2019): A few articles relate to the p-value and alternative
• American Statistician (2019 Supl.): Statistical Inference in the 21st Century: A World Beyond p < 0.05
• Wasserstein (2019): Moving to a World Beyond “p < 0.05”
• McShane (2019): Abandon Statistical Significance
• Amrhein (2019): Comment of Nature: Scientists rise up against statistical significance,
• Editorial of Nature (2019): It’s time to talk about ditching statistical significance
• Different opinions
• Adams (2019): a trillion P values and counting
• Zhang (2019): P values akin to ‘beyond reasonable doubt’
• 2017: Robert Matthews published an article about the changes of statistical practice after ASA's statement on the statistical significance and p-values: The ASA's p-value statement, one year on. I agree on one highlight in this article: "It should be possible to establish firm general principles which focus on what is right rather than what is wrong"
• Abstract:Its aim was to stop the misuse of statistical significance testing. But Robert Matthews argues that little has changed in the 12 months since the ASA's intervention.
• Why do people use p-values instead of computing probability of the model given data?
• 2016: ASA releases statement on statistical significance and p-values (03/07/2016). The statement's six principles:
1. P-values can indicate how incompatible the data are with a specified statistical model.
2. P-values do not measure the probability that the studied hypothesis is true, or the probability that the data were produced by random chance alone.
3. Scientific conclusions and business or policy decisions should not be based only on whether a p-value passes a specific threshold.
4. Proper inference requires full reporting and transparency.
5. A p-value, or statistical significance, does not measure the size of an effect or the importance of a result.
6. By itself, a p-value does not provide a good measure of evidence regarding a model or hypothesis.
• p-value is the likelihood of null hypothesis (H0), a conditional probability given H0:
• p(X ≥ x|H0) for a right tail event
• p(X ≤ x|H0) for a left tail event
• 2 * min{P(X ≥ x|H0), p(X ≤ x|H0)} for a 2-tail event
• α level, the level of significance, is a pre-defined probability of falsely rejecting the null hypothesis we accept. Typically, the statistical significance means the p-value < α level at 0.05.
• False discovery rate (Wikipedia)
• Many times, the p-value was incorrectly interpreted as one of posterior probabilities given the observed data, probabilities given the observed data, the false discovery rate (Colquhoun, 2014) or false positive rate.
• Storey (2003): The Positive False Discovery Rate: A Bayesian Interpretation and the q-Value
• Type I error (α) = 1 - specificity, Type II error (β) = 1 - sensitivity, Power = sensitivity = 1 - β
• Leek (2017). Five ways to fix statistics (Nature)
• Benjamin (2017). Redefine statistical significance
• Baker (2016). Statisticians issue warning over misuse of P values (pdf)
• Matloff (2016). 1) After 150 years, the ASA says no to p-value, 2) Further comments on the ASA manifesto, 3) P-values: the continuing saga
• Benjamini & Galili (2016). It's not the p-values' fault reflections on the recent ASA statement
• Kass (2016). Ten Simple Rules for Effective Statistical Practice. It acknowledged three websites:
• xkcd.com "for conveying statistical ideas with humor"
• Simply Statistics "for thoughtful commentary"
• FiveThirtyEight "for bringing statistics to the world (or at least to the media)".
• Aschwanden (2016). Science Isn’t Broken - It’s just a hell of a lot harder than we give it credit for
• Halsey (2015). The fickle P value generates irreproducible results
• Lazzeroni (2016). Solutions for quantifying P value uncertainty and replication power, and response of Halsey;
• Nuzzo (2013). Scientific method: Statistical errors
• Epimonitor (2016). Growing Concern About Statistical Errors Triggers Statement of P-Values
• Youngquist (2012). Part 19: What is a P value;
• GraphPad's Advice: how to interpret a small P value
• Frost (2014). How to Correctly Interpret P Values, Five Guidelines for Using P values
• Held (2010). A nomogram for P values
• Cumming (2012). Mind your confidence interval: how statistics skew research results
• Gumming (2013). The problem with p values: how significant are they, really?
• Ioannidis (2015). Why most published research findings are false
• Ranstam (2012). Why the P value culture is bad and confidence intervals a better alternative. Ranstam (2009) Sampling uncertainty in medical research. Austin (2002). A brief note on overlapping confidence intervals
• Aschwanden. Statisticians found one thing they can agree on: it's time to stop misusing P-values.
• Lew (2013). Give p a chance: significance testing is misunderstood
• Sullivan (2012). Using Effect Size—or Why the P Value Is Not Enough.
• Fraser (2016). Crisis in Science? or Crisis in Statistics! Mixed messages in Statistics with impact on Science
• Gelman (2014). Data-dependent analysis—a “garden of forking paths”—
• Capital of Statistics. 美国统计协会开始正式吐槽（错用）P值啦
• Wikipedia. Type I & II errors, sensitivity & specificity, effect size, Bayes factor.
• Hubers (2013). Measures of effect size in Stata 13
• Robert Coe (2002). It's the Effect Size, Stupid
• Andrew Gelman (2013). P value and statistical practice, Misunderstanding the p-value
• Simonsohn (2013). Just Post It: The Lesson From Two Cases of Fabricated Data Detected by Statistics Alone
• Goodman (2001). Of P-values and Bayes: a modest proposal
• Goodman (1999). Toward evidence-based medical statistics: the P value fallacy and The Bayes factor (notes: I like the Kass's formula, which uses the likelihood of alternative hypothesis as a numerator, and gives a BF without many decimals).  
• Kass (1995). Bayes factors [on the basis of observed data D, for the dichotomous conditions / models / hypotheses (H1 or H0), Bayes factor = p(D|H1)/p(D|H0)], the rules of thumb assess the quality of the evidence favoring one hypothesis over another as a reference:
• 1 to 3    (not worth more than a bare mention)
• 3 to 20   (positive)
• 20 to 150 (strong)
• > 150     (very strong)
• The odds form of Bayes's theorem for two hypotheses is convenient for calculating a Bayesian update of a chance.
• If there are mutually exclusive hypotheses H1 (= Alternative hypothesis, H1 = Disease) and H0 (= Null hypothesis = No disease) by given D (= Observed data/evidence = Postive test),
• p(H1 and D) = p(H1|D)*p(D) = p(D|H1)*p(H1)
• p(H1|D) = p(H1)*p(D|H1)/p(D)
• p(H0|D) = p(H0)*p(D|H0)/p(D)
• p(H1|D)/p(H0|D) = p(H1)/p(H0)*p(D|H1)/p(D|H0)
• p(H1|D)/[1-p(H1|D)] = p(H1)/[1-p(H1)]*p(D|H1)/p(D|H0)
• odds(H1|D) = odds(H1)*p(D|H1)/p(D|H0)
• Bayes factor (BF) = p(D|H1)/p(D|H0)
• Prior odds of H1 = odds(H1) = p(H1)/p(H0) = p(H1)/[1-p(H1)]
• Posterior odds of H1 given data = odds(H1|D) = p(H1|D)/p(H0|D)= (prior odds)*(BF or likelihood ratio)
• p(H1)=odds(H1)/[1+odds(H1)], p(H1|D)=odds(H1|D)/[1+odds(H1|D)]
• Suppose a person with disease had 3/4 possiblity of positive test, and a person without disease had 1/5 possibility of positive test. and we have no idea of p(D) of that person (the population), which means the prior proability is 50% and the prior odds(D) = p(D) / (1 - p(D) = 1:1. When a person had a positive test:
• the posterior odds (odds(D|Pos)) = (1:1) * (3/4 / 1/5) = 15/4 = 3.75 = (p(D) / (1-p(D)) * (p(Pos|D) / (1 - p(Pos|D)) = odds(D) * BF
• or, the probability had a disease given a positive test (p(D|Pos)) = ((3/4) * 0.5) / (0.5 * 3/4 + (1 - 0.5) * 1/5) = 15/19 = odds(D|Pos) / (1 + odds(D|Pos)) = .79
• Eight versions of Bayes's Theorem (pdf): simple, explicity, general, Sigma, canceled, odds, relative odds, and compound odds.
• Currell (2009). Chapter 7 Bayesian statistics
• O’Hagan (2006). Bayes factors, Deeks (2004). Diagnostic tests 4: likelihood ratios, Lindley (2004). Bayesian thoughts, Griffis (2006). Statistics and the Bayesian mind
• Berger (1988). The likelihood principle: a review, generalizations, and statistical implications
• Masson (2011): A tutorial on a practical Bayesian alternative to null-hypothesis significance testing 
• Faulkenberry (2018): A Simple Method for Teaching Bayesian Hypothesis Testing in the Brain and Behavioral Sciences
• Sellke (2001): Calibration of p Values for Testing Precise Null Hypotheses
• Berger (1987): Testing a Point Null Hypothesis: The Irreconcilability of P Values and Evidence
• Bayesian p-value
• Nowozin (2015). Bayesian P-Values
• Lin (2009). Using Bayesian p-values in a 2 x 2 table of matched pairs with incompletely classified data
• NISS Webinar (2019): Alternatives to the Traditional P-value
  

Curve-Fitting

XKCD: Curve-Fitting

how to estimate the risk and relative risk from logistic regression

How to estimate the risk and relative risk from logistic regression of a case-control study (OR to RR)
One advantage of odds ratios is that we can estimate it in a case-control study, when we usually oversample the cases. However, we cannot directly calculate the probability using logistic regression in a case-control study, since the beta(0) of the case-control study cannot represent the target population, which beta(0) of the target population is equal to beta(0) of the case-control study - log(sampling probability for cases/sampling probability for controls).

• King (2004). Inference under full population information in "Case-Control Studies, Inference in" of Encyclopedia of Biopharmaceutical Statistics
• Zhang (1998). What's the relative risk? A method of correcting the odds ratio in cohort studies of common outcomes
• Greenland (1981). Multivariate estimation of exposure-specific incidence from case-control studies
• Austin (2010). Overestimation of risk ratios by odds ratios in trials and cohort studies: alternatives to logistic regression
• Austin (2011). A Tutorial on Methods to Estimating Clinically and Policy-Meaningful Measures of Treatment Effects in Prospective Observational Studies: A Review
• Spiegelman (2005). Easy SAS calculations for risk or prevalence ratios and differences
• Bartlett (2015). Estimating risk ratios from observational data in Stata
• Chui (2013). A Regression-Based Method for Estimating Risks and Relative Risks in Case-Base Studies
• Knol (2012). Overestimation of risk ratios by odds ratios in trials and cohort studies: alternatives to logistic regression

R functions and keyboard shortcuts

R functions/commands and keyboard shortcuts

R! is powerful and has rich packages and functions. It's impossible to build a list of functions/shortcuts to fit the purposes of all. Below are some functions/shortcuts related to my projects.

• Cheatsheets, The R Guide, R Reference Card
• Help functions: help()/?, apropos(), find(): apropos() finds all objects. find() the locations of found objects, methods(), example(), demo(), vignette(), args() 
• Housekeeping functions: getwd(), setwd(), rm(list=ls()) removes all objects in the R environment, source("myRscript.r") runs the R codes in "myRscript.r" file, fix() modifies the original object, and edit() is used edit an object and returns to a new object, download.file() downloads a file from the Internet, attach()/detach() objects, search() shows the current search paths and sequence, install.packages(), update.packages(), remove.packages(), getOption("defaultPackages") which can be changed by setting the option in startup code (e.g. in ~/.Rprofile), .libPaths()
• Numeric/character functions: length(), seq(), rep(), cut(), pretty(), cat(), substr(), grep(), sub(), strsplit(), paste(), toupper(), tolower()
• Data functions: read.table(), head(), tail(), str(), class(), length(), dim(), nrow(), ncol(), names(), levels(), length(), c(), cbind(), rbind(), append(), rep(), rev(), sort(), unique()
• Type functions: "is." for checking or "as." for conversion + numeric(), character(), vector(), matrix(), data.frame(), factor(), logical(), integer(). For example: is.numeric(), as.numeric()
• Mathematical functions: abs(), sqrt(), log(), log(x, base=n), log10(), exp(), prod(), factorial(), choose(), ceiling(), floor(), solve(), trunc(), round(), signif(), cos(), sin(), tan(), acos()
• Statistical functions: mean(), median(), sd(), var(), mad(), quantile(), range(), sum(), diff(), min(), max(), scale(), fivenum(), cumsum(), cumprod(), cummax(), cumin(), cor(), colSums(), rowSums(), colMeans(), rowMeans()
• Probability functions: the form is [d][p][q][r]distribution(). d, p, q, r are for (d)ensity, cumulated (p)robability/distribution function, (q)uantile function, and (r)andom generation, respectively. the Distribution types can be: (norm)al, (beta), (binom)ial, (chisq)uared, (exp)onential, (logis)tic, (multinom)ial, (n)egative (binom)ial, (pois)son, (f), (gamma), (t), (unif)orm, etc. for example: dnorm(), pnorm(), qnorm(), rnorm()
• Statistical modeling functions
• Model functions: lm(), glm(), nls(), nls2(), lme() / nlme()
• Symbol formulas (y ~ A + B + C ): ":" is for interaction term, "*" is for complete interaction, "^" is for crossing to a specified degree "." is a placeholder for all other variables except the dependent variable, "-" removes a variable from the equation, "-1" suppresses the intercept, "I()" has elements within the parentheses interpreted arithmetically
• Post-estimation functions: coef(), confint(), resid(), fitted(), summary(), predict(), deviance(), print(),plot(), formula(), anova(obj1, obj2), AIC(), vcov()
• Contrast functions: contr.helmert(), contr.poly(), contr.sum(), contr.treatment(), contr.SAS()
• RStudio is an integrated development environment (IDE) for R. RStudio combines an intuitive user interface with powerful coding tools to help you get the most out of R. Shortcuts (you can modify them: Tools -> Modify Keyboard Shortcuts...)
• Alt + Shift + K: Show a Quick Reference
• Alt + -: Insert assignment operator "<- font="">
• Ctrl + Shift + M: Insert pipe operator "%>%" (I changed it as Ctrl + Shift + P)
• Ctrl + Alt + I: Insert chunk (R Notebook/Markdown)
• Ctrl + 1: Move cursor to source Editor window
• Ctrl + 2: Move cursor to Command window
• Ctrl + 3: Move cursor to Help window
• Ctrl + 4: Move cursor to History window
• Ctrl + 5: Move cursor to File window
• Ctrl + 6: Move cursor to Plots window
• ...

choice of analytical language

Choice of analytical language
I have used mainly three statistical languages, Stata, R, and SAS, for many years for different purposes. The weights of usage of those three languages are shift from SAS-Stata-R to SAS-R-Stata, then, to Stata-R-SAS. Sometimes I am asked to recommend a better analytic language, which is always a hard and complicated question to me. I came across an blog written by Curtis Miller, which is very thoughtful and helpful to make this kind of choice. Here is his blog: "On Programming Languages; Why My Dad Went From Programming to Driving a Bus". Hopefully his story can help you to make your own decision.

Using multi-year national survey cohorts for period estimates: an application of weighted discrete Poisson regression for assessing annual national mortality in US adults with and without diabetes, 2000-2006.

Cheng YJ, Gregg EW, Rolka DB, Thompson TJ.

BACKGROUND:

Monitoring national mortality among persons with a disease is important to guide and evaluate progress in disease control and prevention. However, a method to estimate nationally representative annual mortality among persons with and without diabetes in the United States does not currently exist. The aim of this study is to demonstrate use of weighted discrete Poisson regression on national survey mortality follow-up data to estimate annual mortality rates among adults with diabetes.

METHODS:

To estimate mortality among US adults with diabetes, we applied a weighted discrete time-to-event Poisson regression approach with post-stratification adjustment to national survey data. Adult participants aged 18 or older with and without diabetes in the National Health Interview Survey 1997-2004 were followed up through 2006 for mortality status. We estimated mortality among all US adults, and by self-reported diabetes status at baseline. The time-varying covariates used were age and calendar year. Mortality among all US adults was validated using direct estimates from the National Vital Statistics System (NVSS).

RESULTS:

Using our approach, annual all-cause mortality among all US adults ranged from 8.8 deaths per 1,000 person-years (95% confidence interval [CI]: 8.0, 9.6) in year 2000 to 7.9 (95% CI: 7.6, 8.3) in year 2006. By comparison, the NVSS estimates ranged from 8.6 to 7.9 (correlation = 0.94). All-cause mortality among persons with diabetes decreased from 35.7 (95% CI: 28.4, 42.9) in 2000 to 31.8 (95% CI: 28.5, 35.1) in 2006. After adjusting for age, sex, and race/ethnicity, persons with diabetes had 2.1 (95% CI: 2.01, 2.26) times the risk of death of those without diabetes.

CONCLUSION:

Period-specific national mortality can be estimated for people with and without a chronic condition using national surveys with mortality follow-up and a discrete time-to-event Poisson regression approach with post-stratification adjustment. (Full text)

accept-reject algorithm

Accept-reject algorithm
Accept-reject algorithm (acceptance-rejection method) or reject sampling is a simple and general simulation method to decide observations with or without a trait from the probability of a distribution. In this way, we can convert a probability into a dichotomous condition (i.e. yes or no). Basically, there are three steps:

• Step 1. Generate Y from density g [Y = f(x), the pdf of f(x) is the target distribution]
• Sample a point (an x-position) from the proposal density distribution (g) and draw a vertical line at this point, get the density (an y-position) [X ~ g(x)]. The density function of Y has a upper, a constant c, and c is >=1.
• Step 2. Generate U from the uniform distribution on the interval (0, cg(x)) [U = cg(x), the pdf of cg(x) is the proposal distribution]
• Sample uniformly along in the range of x-position (i.e. uniformly from 0 to the maximum of the probability density function) [U ~ runif(0, 1)]
• Step 3. If U <= Y, then set Y = X ("accept"), else repeat Steps 1 and 2

Pr(X|accept) = Pr(accept|X) x Pr(X)/Pr(accept), using Bayes' theorem
Pr(accept|X) = f(x)/cg(x)
Pr(X) = g(x)
Pr(accept) = 1/c
therefore, Pr(X|accept) = f(x)


Example: Stata simulation and define the event

clear
set seed 770488
set obs 1000

gen x = runiform() - .5
gen z = runiform() - .5
gen xb = x + 8*z
 gen y = 1 / (1 + exp(xb)) < runiform() // y defined as 0 or 1
logistic y x z

general linear models vs. generalized linear models

General linear models vs. generalized linear models

• Wikipedia: Comparison of general and generalized linear models (modified)

	General linear model	Generalized linear model
Typical estimation method	Least squares, best linear unbiased prediction	Maximum likelihood or Bayesian
Special cases	ANOVA, ANCOVA, MANOVA, MANCOVA, ordinary linear regression, mixed model, t-test, F-test	linear regression, logistic regression, Poisson regression, gamma regression
Function in R	lm()	glm()
Function in Matlab	mvregress()	glmfit()
Procedure in SAS	PROC GLM, PROC MIXED	PROC GENMOD (PROC LOGISTIC for logistic regression only), PROC GLIMMIX Comparing the MIXED and GLIMMIX
Command in Stata	regress	Stata Generalized linear models
Function in Mathematica	LinearModelFit	GeneralizedLinearModelFit
Command in EViews	ls

• Generalized linear models have the flexiblility for response variables that have other than a normal distribution. If a generalized linear model uses an identity link function and a normal family distribution, then this model is equivalent to a general linear model.
• Generalized linear mixed models have the flexibility to model random effects and correlated errors for nonmormal data.