MIDAS Regression

[Vishv Jeet]

I am trying to set up a univariate linear regression in Stan. Consider the following Stan program:

data {

  int N;

  int K;

  vector[N] y;

  matrix[N, K] X;

}

parameters {

  real<lower = 0> alpha;

  real<lower = 0> beta;

  real<lower = 0> lambda;

  real<lower = 0> sigma;

}

transformed parameters{

  vector<lower = 0>[K] weights;

  vector[N] wX;

  for (k in 1:K) {

    weights[k] = exp(beta_lpdf((k + 0.0) / K| alpha, beta));

  }

  wX = X * weights / sum(weights);

}

model {

  y ~ normal(wX * lambda, sigma);

  alpha ~ cauchy(0, 2);

  beta ~ cauchy(0, 2);

  lambda ~ cauchy(0, 2);

  sigma ~ cauchy(0, 2);

}

Now in addition to lambda and sigma, the parameters that an OLS style regression between y and wX can estimate, I also want to estimate alpha and beta of a Beta distribution of weights with which I reduce matrix X to a single vector wX.

I think I have set it up correctly except the use of beta_lpdf function. I am getting the following error from Stan:

[1] "Rejecting initial value:"                                                                             

  [2] "  Gradient evaluated at the initial value is not finite."                                             

  [3] "  Stan can't start sampling from this initial value."   

Any help/comments/suggestions are much appreciated.

[Vishv Jeet]

The following program works, the B(alpha, beta) cancels out when normalizing. I think the call to beta_lpdf was wrong but I would still like to understand that better.

data {

  int N;

  int K;

  vector[N] y;

  matrix[N, K] X;

}

parameters {

  real<lower = 0> alpha;

  real<lower = 0> beta;

  real<lower = 0> lambda;

  real<lower = 0> sigma;

}

transformed parameters{

  vector<lower = 0>[K] weights;

  vector[N] wX;

  {

     real x;

     for (k in 1:K) {

        x = k + 0.0) / K;

        weights[k] = x^alpha * (1 - x)^beta;

     }

[Bob Carpenter]

In the second program you have:

> for (k in 1:K) {
> x = k + 0.0) / K;
> weights[k] = x^alpha * (1 - x)^beta;
> }

Which can't possibly compile due to the missing open paren defining x.

In the first you have:

> for (k in 1:K) {

> weights[k] = exp(beta_lpdf((k + 0.0) / K| alpha, beta));
> }

Those don't give you the same answer even if you fix the syntax.
If you look at the definition of the beta_lpdf you'll see that
there are gamma terms to normalize which aren't included in
your kernel-style version in the second program.

I don't know what you're trying to do---just saying they're not
the same program.

You can test this by just writing the function out and defining
the weights as generated quantities (or transformed parameters). You'd
see they are different.

- Bob

> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

[John Hall]

It looks to me like you're putting a sampling statement into a variable, which doesn't make much sense.

I think what you would have wanted (based on https://pdfs.semanticscholar.org/4d02/f73c1acf1af08090639f166ccde8f691e2bc.pdf) was something like
weights[k] = x^alpha * (1 - x)^beta / exp(lbeta(alpha, beta));
which might be equivalent to
weights[k] = inc_beta(alpha, beta, x);
but I would defer to someone else to confirm.

Anyway, I think the first one makes it pretty obvious why it cancels out when you divide by the sum.

You could also speed things up by doing the calculation of x in transformed data. You might also be able to vectorize the weights calculation.