Table of Contents

Preface


1. Astrostatistics

"Begin at the beginning," the King said, very gravely, "and go on till you come to the end: then stop". 

Alice in Wonderland,  Lewis Carroll (1865)


2. Prerequisites

             Don’t Panic.
             The Hitchhiker’s Guide to the Galaxy, Douglas Adams (1981)

 

Code 2.1  - Example of linear regression in R

Code 2.3  - Example of linear regression in Python

 


3. Frequentist vs Bayesian methods

It is a capital mistake to theorize before one has data. Insensibly one  begins to twist facts to suit theories, instead of theories to suit facts.
                                                                   Sherlock Holmes, Arthur Conan Doyle (1887)

Code 3.1  - Basic linear model in R

Code 3.2  - Ordinary least squares regression in Python without formula

Code 3.3  - Bayesian normal linear model in R

Code 3.4  - Bayesian normal linear model in Python


4. Normal linear models

Let us then suppose the mind to be, as we say, white paper, void of all characters, without any ideas: – How comes it to be furnished? Whence comes it by that vast store which the busy and boundless fancy of man has painted on it with an almost endless variety? To this I answer, in one word, from EXPERIENCE.
                                                                                                     
 John Locke (1698)

Code 4.1  - Normal linear model in R using JAGS

Code 4.2  - Normal linear model in R using JAGS and the zero trick

Code 4.3  - Normal linear model in Python using Stan

Code 4.4  - Plotting posteriors from Code 4.3

Code 4.6  - Multivariate normal linear model in R using JAGS

Code 4.7  - Multivariate normal linear model in Python using Stan

Code 4.8 -  Synthetic normal data in R with errors in measurements

Code 4.9  - Normal linear model in R using JAGS and ignoring errors in measurements

Code 4.10  - Normal linear model in R using JAGS and including errors in variables

Code 4.11  - Normal linear model in Python using Stan and including errors in variables


5. GLM part I - continuous and binomial models

Sit down before a fact as a little child, be prepared to give up every preconceived notion. Follow humbly wherever and to whatever abysses nature leads, or you shall learn nothing.
                                                                                            
Thomas H. Huxley (1860)

Code 5.1 - GLM logistic regression in R

Code 5.2 - GLM logistic regression in Python

Code 5.3 - Synthetic lognormal data generated in R

Code 5.4 - Lognormal model in R using JAGS

Code 5.5 - Lognormal model in Python using Stan

Code 5.6 - Log-gamma synthetic data generated in R

Code 5.7 - Log-gamma model in R using JAGS

Code 5.8 - Log-gamma model in Python using Stan

Code 5.9 - Log-inverse-Gaussian data

Code 5.10 - Log-inverse-Gaussian model in R using JAGS

Code 5.11 - Inverse Gaussian model in Python using Stan

Code 5.12 - Synthetic beta-distributed data generated in R

Code 5.13 - Beta model in R using JAGS

Code 5.14 - Beta model in Python using Stan

Code 5.15 - Synthetic data from logistic model in R

Code 5.16 - Logistic model using R

Code 5.17 - Logistic model in R using JAGS

Code 5.18 - Logistic model using pymc3

Code 5.19 - Logistic model in Python using Stan

Code 5.20  -Synthetic probit data and model generated in R

Code 5.21 - Probit model using R

Code 5.22 - Probit model in R using JAGS

Code 5.23 - Probit model in Python

Code 5.24 - Probit model in Python using Stan

Code 5.25 - Synthetic data from a binomial model in R

Code 5.26 - Binomial model in R using JAGS

Code 5.28 - Binomial model in Python using Stan

Code 5.31 - Simulated beta–binomial data in R

Code 5.32 - Beta–binomial synthetic model in R using JAGS

Code 5.33 - Explicitly given beta–binomial data in R

Code 5.34 - Beta–binomial model in R using JAGS for explicitly given data and the zero trick

Code 5.35 - Beta–binomial model with inverse link in R using JAGS

Code 5.36 - Beta–binomial model with synthetic data in Python using Stan


6. GLM part II - count models

And now, each night I count the stars.
And each night I get the same number.
And when they will not come to be counted,
I count the holes they leave.

Preface to a Twenty Volume Suicide Note, Amiri Baraka (1961)

Code 6.1 - Synthetic data following a Poisson distribution in R

Code 6.2 - Synthetic Poisson data and model in R: binary and continuous predictors

Code 6.3 - Bayesian Poisson model using R

Code 6.4 - Bayesian Poisson model in R using JAGS

Code 6.5 - Simple Poisson model in Python

Code 6.6 - Multivariate Poisson model in Python

Code 6.7 - Bayesian Poisson model using pymc3

Code 6.8 - Bayesian Poisson model in Python using Stan

Code 6.9 - Synthetic negative binomial data and model in R

Code 6.10 - Negative binomial model in R using COUNT

Code 6.11 - Bayesian negative binomial in R using JAGS

Code 6.13 - Negative binomial: direct parameterization using JAGS and dnegbin

Code 6.14 - Negative binomial with zero trick using JAGS directly

Code 6.16 - Negative binomial model in Python using pymc3

Code 6.17 - Negative binomial model in Python using Stan

Code 6.18 - Synthetic data for generalized Poisson

Code 6.19 - Bayesian generalized Poisson using JAGS

Code 6.20 - Generalized Poisson model in Python using Stan

Code 6.21 - Zero-truncated Poisson data

Code 6.22 - Zero-truncated Poisson with zero trick

Code 6.23 - Zero-truncated Poisson model in Python using Stan

Code 6.24 - Zero Truncated Negative binomial with 0-trick using JAGS - direct

Code 6.25 - Zero-truncated negative binomial model in Python using Stan

Code 6.26 - Create synthetic negative binomial data

Code 6.27 - Bayesian three-parameter NB-P – indirect parameterization with zero trick

Code 6.28 - Negative binomial model with three parameters in Python using Stan. 


7. GLM part III - zero-inflated and hurdle models

 Nothing isn’t better or worse than anything. Nothing is just nothing.
                           
Arya Stark - A Game of Thrones,  George R. R. Martin (1996)

Code 7.1 - Bayesian zero-inflated Poisson model in R using JAGS

Code 7.2 - Bayesian zero-inflated Poisson in Python using Stan

Code 7.3 - Zero-inflated negative binomial synthetic data in R

Code 7.4 - Bayesian zero-inflated negative binomial model using JAGS

Code 7.5 - Bayesian zero-inflated negative binomial model in Python using Stan

Code 7.6 - Synthetic data for Poisson–logit hurdle zero-altered models

Code 7.7 - Bayesian Poisson–logit hurdle zero-altered models

Code 7.8 - Bayesian Poisson–logit hurdle model in Python using Stan

Code 7.9 - Zero-altered negative binomial (ZANB) or NB hurdle model in R using JAGS

Code 7.10 - Zero-altered negative binomial (ZANB) or NB hurdle model in Python using Stan

Code 7.11 - Bayesian log-gamma–logit hurdle model in R using JAGS

Code 7.12 - Bayesian log-gamma–logit hurdle model in Python using Stan

Code 7.13 - Bayesian lognormal–logit hurdle using JAGS

Code 7.14 - Bayesian lognormal–logit hurdle model in Python using Stan


8. Hierarchical GLMMs

       The complexity of things, the things within things - just seems to be endless.
                                                                                                   
 Alice Munro (2001)

Code 8.1 - Random intercept Gaussian data generated in R

Code 8.2 - Random intercept normal model in R using JAGS

Code 8.4 - Random intercept normal model in Python using Stan

Code 8.5 - Simulated random intercept binary logistic data

Code 8.6 - Bayesian random intercept binary model in R

Code 8.7 - Bayesian random intercept binary logistic model in Python using pymc3

Code 8.8 - Bayesian random intercept binary logistic model in R using JAGS

Code 8.9 - Random intercept binomial logistic data in R

Code 8.10 - Random intercept binomial logistic model in R using JAGS

Code 8.11 - Random intercept binomial logistic model in Python using Stan

Code 8.13 - Bayesian random intercept Poisson model in Python

Code 8.14 - Bayesian random intercept Poisson model in R using JAGS

Code 8.15 - Bayesian random intercept Poisson model in Python using Stan

Code 8.16 - Random-intercept–random-slopes Poisson data in R

Code 8.17 - Random-intercept–random-slopes Poisson model in R using JAGS

Code 8.18 - Random-intercept–random-slopes Poisson model in Python using Stan

Code 8.19 - Random Intercept negative binomial data in R

Code 8.20 - Random intercept negative binomial mixed model in R

Code 8.21 - Bayesian random intercept negative binomial mixed model in Python using pymc3

Code 8.22 - Bayesian random intercept negative binomial in R using JAGS

Code 8.23 - Bayesian random intercept negative binomial in Python using Stan


9. Model selection

Crying is all right in its way while it lasts.

But you have to stop sooner or later, and then you still have to decide what to do.

 C.S. Lewis - The Silver Chair (Chronicles of Narnia, #4, 1953)

Code 9.4 - K and M model in Python using Stan


10. Astronomical applications

       MAKE IT SO.
       Captain Jean-Luc Picard,  Star Trek: The Next Generation (1987)

 

 

Code 10.1 - Normal linear model in R using JAGS for accessing the relationship between central black hole

                    mass and bulge velocity dispersion

Code 10.2 - Normal linear model, in Python using Stan, for assessing the relationship between central black

                    hole mass and bulge velocity dispersion

Code 10.3 - Gaussian linear mixed model, in R using JAGS, for modeling the relationship    between type Ia                      supernovae host galaxy mass and Hubble residuals

Code 10.4 - Gaussian linear mixed model, in Python using Stan, for modeling the relationship between type

                    Ia supernovae host galaxy mass and Hubble residuals

Code 10.5 - Multivariate    normal   model    in    R using    JAGS   for    accessing the relationship between

                    period, luminosity, and color in early-type contact binaries

Code 10.6 - Multivariate Gaussian mixed model in Python, using Stan, for accessing the relationship

                    between  luminosity, period, and color in early-type contact binaries

Code 10.7 - Lognormal model in R using JAGS to describe the initial mass function (IMF)

Code 10.8 - Plotting routine, in R, for Figure 10.6

Code 10.9 - Lognormal model in Python using Stan to describe the initial mass function (IMF)

Code 10.10 - Beta model in R using JAGS, for accessing the relationship   between   the   baryon   fraction 

                      in atomic gas and galaxy stellar mass

Code 10.11 - Beta model in Python using Stan, for accessing the  relationship  between the fraction of 

                      atomic gas and the galaxy stellar mass

Code 10.12 - Bernoulli model in R  using JAGS,  for  accessing  the  relationship  between  bulge size and

                      the fraction of red spirals

Code 10.13 - Bernoulli model in Python using Stan, for assessing the relationship between bulge size and 

                      the fraction of red spirals.

Code 10.14 - Poisson model, in R using JAGS, for modeling the relation between globular clusters

                      population and host galaxy visual magnitude

Code 10.15 - Negative binomial model, in R using JAGS, for   modeling  the  relationship  between  

                      globular cluster population and host galaxy visual magnitude

Code 10.16 - NB-P   model   in   R   using JAGS,   for   modeling   the  relationship  between  globular 

                      cluster population and host galaxy visual magnitude

Code 10.17 - Negative binomial model in Python using Stan, for modeling the relationship between globular

                      cluster population and host galaxy visual magnitude

Code 10.18 - Bernoulli logit model, in R using JAGS, for accessing   the   relationship  between Seyfert

                      AGN activity and galactocentric distance

Code 10.19 - Bernoulli logit model, in Python using Stan, for assessing the relationship between Seyfert

                      AGN activity and galactocentric distance

Code 10.20 - Lognormal–logit hurdle model, in R using JAGS, for assessing  the  relationship  between

                      dark-halo mass and stellar mass

Code 10.21 - Lognormal–logit hurdle model, in Python  using  Stan,  for  assessing  the  relationship

                      between dark halo mass and stellar mass

Code 10.22 - Normal autoregressive  model  AR(1) for  accessing the  evolution  of  the  number  of  

                      sunspots through the years.

Code 10.23 - Model shown in Figure 10.21

Code 10.24 - Negative binomial model (AR1) for assessing the evolution of the  number  of  sunspots

                      through the years.

Code 10.25 - Negative binomial model (AR1) in Python using Stan, for assessing the evolution of the

                      number of sunspots through the years

Code 10.26 - Bayesian normal model for cosmological parameter inference from type Ia supernova  data in

                      R  using Stan.


11. The future of Astrostatistics

There is nothing like a dream to create the future.
                                    Les Miserables, Victor Hugo (1862)

Appendix A. Bayesian modeling using INLA

Bibliography

Index

© 2017 by Emille E. O. Ishida