Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Document binomial-logit GLM #678

Merged
merged 3 commits into from
Jan 16, 2024
Merged
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
126 changes: 126 additions & 0 deletions src/functions-reference/bounded_discrete_distributions.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -155,6 +155,132 @@ The log binomial probability mass of n successes in N trials given
logit-scaled chance of success alpha dropping constant additive terms
`r since("2.25")`

## Binomial-logit generalized linear model (Logistic Regression) {#binomial-logit-glm}

Stan also supplies a single function for a generalized linear model
with Binomial likelihood and logit link function, i.e. a function
andrjohns marked this conversation as resolved.
Show resolved Hide resolved
andrjohns marked this conversation as resolved.
Show resolved Hide resolved
for logistic regression with aggregated outcomes. This provides a more efficient
implementation of logistic regression than a manually written
regression in terms of a Binomial likelihood and matrix
multiplication.

### Probability mass function

Suppose $N \in \mathbb{N}$, $x\in \mathbb{R}^{n\cdot m}, \alpha \in \mathbb{R}^n, \beta \in \mathbb{R}^m$, and $n \in
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would start by saying that M is the number of predictors and N is the number of data items.

Then we need the notation to match with caps and to use \times not \cdot

  • x is in R^(M x N) (with \times, not \cdot)
  • \alpha in \mathbb{R}
  • \beta \in \mathbb{R}^M

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Then we need the notation to match with caps and to use \times not \cdot

That's also inconsistent with the other _glm docs:

\{0,\ldots,N\}$. Then \begin{align*}
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Start the begin{align*} on its own line and line up the text under it so that it's readable. Looks like this accidentally got collapsed into a paragraph.

In the Stan doc, if we have an M x N matrix, we've conventionally used [m, n] for indexing, not [i, j].

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In the Stan doc, if we have an M x N matrix, we've conventionally used [m, n] for indexing, not [i, j].

That's not the case for the _glm docs

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Then we should change the _glm docs to be consistent with the rest of the User's Guide. I never wrote down hard and fast style rules and things like the symbol for expectation starts drifting under multiple authors. I mean to go and do a consistency fix of the entire doc set soon.

&\text{BinomialLogitGLM}(n~|~N, x, \alpha, \beta) = \text{Binomial}(n~|~N,\text{logit}^{-1}(\alpha_i + x_i\cdot
\beta)) \\ &= \binom{N}{n} \left( \text{logit}^{-1}(\alpha_i + \sum_{1\leq j\leq m}x_{ij}\cdot \beta_j) \right)^{n} \left( 1 -
\text{logit}^{-1}(\alpha_i + \sum_{1\leq j\leq m}x_{ij}\cdot \beta_j) \right)^{N - n}. \end{align*}

### Sampling statement

`n ~ ` **`binomial_logit_glm`**`(N, x, alpha, beta)`

Increment target log probability density with `binomial_logit_glm_lupmf(n | N, x, alpha, beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm ~; -->
\index{{\tt \bfseries binomial\_logit\_glm }!sampling statement|hyperpage}

### Stan Functions

<!-- real; binomial_logit_glm_lpmf; (int n | int N, matrix x, real alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lpmf }!{\tt (int n \textbar\ int N, matrix x, real alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lpmf`**`(int n | int N, matrix x, real alpha, vector beta)`<br>\newline
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Like our other function doc, this should name the arguments. alpha is an intercept, beta is a vector slopes, x is the data matrix, N is the total count and n is the count of successes.

Now if x is a matrix, doesn't n and N have to be 1D arrays?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Like our other function doc, this should name the arguments. alpha is an intercept, beta is a vector slopes, x is the data matrix, N is the total count and n is the count of successes.

That would be inconsistent with the current doc for other _glm functions:

Now if x is a matrix, doesn't n and N have to be 1D arrays?

No, they would be broadcast to match. This is the same way in the bernoulli_logit_glm doc

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Then normal_id_glm and bernoulli_logit_glm should be fixed to match the rest of our doc. Given that we've thrown out consistency, there's no need for any new GLM to be consistent with the other GLM doc. I'd rather it be consistent with the rest of our doc as that's where it's going to go in the end.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I should've added that this doesn't need to be done as part of this PR. If you leave the new _glm like the other GLM code, I can just fix it in a pass to make everything consistent again.

The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm_lupmf; (int n | int N, matrix x, real alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lupmf }!{\tt (int n \textbar\ int N, matrix x, real alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lupmf`**`(int n | int N, matrix x, real alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)` dropping constant additive terms.
`r since("2.34")`

<!-- real; binomial_logit_glm_lpmf; (int n | int N, matrix x, vector alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lpmf }!{\tt (int n \textbar\ int N, matrix x, vector alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lpmf`**`(int n | int N, matrix x, vector alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Don't need to say "mass" here---it's just a probability. Or it's a "log probability mass function" if you want to spell it all out.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

"probability mass" is used throughout the other function docs:

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'd argue that "probability mass" is not idiomatic in this context because we just call the resulting quantity "probability" not a "probability mass". It's a "probability mass function" but it returns a probability not a probability mass. So this is another case where the binomial, Bernoulli, etc. need to be fixed. No worries if you can't get to it in this PR.

`inv_logit(alpha + x * beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm_lupmf; (int n | int N, matrix x, vector alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lupmf }!{\tt (int n \textbar\ int N, matrix x, vector alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lupmf`**`(int n | int N, matrix x, vector alpha, vector beta)`<br>\newline
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It doesn't make sense for x to be a matrix and n and N to be scalars. Is there interpretation here that you are broadcasting the n and N for all of the rows of x?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's right, this is consistent with the signatures for bernoulli_logit_glm and normal_id_glm

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks. I think it'd help clarify this in the doc. For instance, line 315 says "The log normal probability density of y given location alpha + x * beta and scale sigma." but in this case y is a vector and alpha + beta * x is a vector, so calling a vector a location seems to violate agreement (plural/singular).

The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)` dropping constant additive terms.
`r since("2.34")`

<!-- real; binomial_logit_glm_lpmf; (array[] int n | array[] int N, row_vector x, real alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lpmf }!{\tt (array[] int n \textbar\ array[] int N, row\_vector x, real alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lpmf`**`(array[] int n | array[] int N, row_vector x, real alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm_lupmf; (array[] int n | array[] int N, row_vector x, real alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lupmf }!{\tt (array[] int n \textbar\ array[] int N, row\_vector x, real alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lupmf`**`(array[] int n | array[] int N, row_vector x, real alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)` dropping constant additive terms.
`r since("2.34")`

<!-- real; binomial_logit_glm_lpmf; (array[] int n | array[] int N, row_vector x, vector alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lpmf }!{\tt (array[] int n \textbar\ array[] int N, row\_vector x, vector alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lpmf`**`(array[] int n | array[] int N, row_vector x, vector alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm_lupmf; (array[] int n | array[] int N, row_vector x, vector alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lupmf }!{\tt (array[] int n \textbar\ array[] int N, row\_vector x, vector alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lupmf`**`(array[] int n | array[] int N, row_vector x, vector alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)` dropping constant additive terms.
`r since("2.34")`


<!-- real; binomial_logit_glm_lpmf; (array[] int n | array[] int N, matrix x, real alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lpmf }!{\tt (array[] int n \textbar\ array[] int N, matrix x, real alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lpmf`**`(array[] int n | array[] int N, matrix x, real alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm_lupmf; (array[] int n | array[] int N, matrix x, real alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lupmf }!{\tt (array[] int n \textbar\ array[] int N, matrix x, real alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lupmf`**`(array[] int n | array[] int N, matrix x, real alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)` dropping constant additive terms.
`r since("2.34")`

<!-- real; binomial_logit_glm_lpmf; (array[] int n | array[] int N, matrix x, vector alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lpmf }!{\tt (array[] int n \textbar\ array[] int N, matrix x, vector alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lpmf`**`(array[] int n | array[] int N, matrix x, vector alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)`.
`r since("2.34")`

<!-- real; binomial_logit_glm_lupmf; (array[] int n | array[] int N, matrix x, vector alpha, vector beta); -->
\index{{\tt \bfseries binomial\_logit\_glm\_lupmf }!{\tt (array[] int n \textbar\ array[] int N, matrix x, vector alpha, vector beta): real}|hyperpage}

`real` **`binomial_logit_glm_lupmf`**`(array[] int n | array[] int N, matrix x, vector alpha, vector beta)`<br>\newline
The log Binomial probability mass of n given N trials and chance of success
`inv_logit(alpha + x * beta)` dropping constant additive terms.
`r since("2.34")`

## Beta-binomial distribution

### Probability mass function
Expand Down