linear regression: infer points (mean, sd) instead of slopes

[Stefano]

Hello Community,

I can't get around this issue.

normally we apply linear regression in this way (with intercept = 0):

vector y
vector x 1,2,3
real beta 1,2,3

y ~ normal(beta1 * x1 + beta2 * x2 + beta3 * x3, sigma);

now, let's suppose alpha and beta1,2,3 are known, for many observations (random numbers for this example)

vector y
real x 1,2,3
vector beta 1,2,3

y	beta1	beta2	beta3
25.87397	0.162388	0.762557	0.075055
43.63101	0.221546	0.339176	0.439278
72.96114	0.187963	0.434091	0.377947
3.900392	0.291026	0.2098	0.499174
23.68992	0.321648	0.137923	0.540429
82.86007	0.144561	0.321849	0.53359
92.71681	0.753214	0.130166	0.116621
72.43547	0.2165	0.195072	0.588428
39.80495	0.473216	0.497003	0.029782
70.54467	0.297066	0.274241	0.428693
50.53164	0.301017	0.310204	0.388779
61.76674	0.136221	0.522591	0.341188
80.84892	0.172414	0.443751	0.383835

How can I infer x1,3,4 both mean and standard deviation? predicting the point values of x1,2,3 it's easy, but what if I want to predict also their standard deviation?

Thanks a lot.

[Bob Carpenter]

You need a model for the missingness. See the chapter on missing data
in the Stan manual. In a typical Bayesian setting, you create a joint
model p(y, x, beta) and then turn the usual Bayesian crank to draw
from p(x | y, beta) just as you would draw from p(beta | x, y).
That is, declare y and beta as data and x as a parameter, code the
joint density in Stan and push the button.

- 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.

[Stefano]

Thanks Bob,

On Thursday, November 10, 2016 at 2:57:48 AM UTC+11, Bob Carpenter wrote:

You need a model for the missingness.  See the chapter on missing data
in the Stan manual.  In a typical Bayesian setting, you create a joint
model p(y, x, beta) and then turn the usual Bayesian crank to draw
from p(x | y, beta)

so you think in principle I can expand doing p(x_mu, x_sd | y, beta)?

Because you proposal works, but I didn't manage to infer SD for each of x 1, 2, 3 (the model does not converge)

 

just as you would draw from p(beta | x, y).  
That is, declare y and beta as data and x as a parameter, code the
joint density in Stan and push the button.

- Bob

 

Thanks

[Bob Carpenter]

Did you put a proper prior on x? What didn't converge?

Yes, you can do any inferences you want, including
p(mu, sigma | y, beta) or whatever, but you'll need priors
on the variables and you may run into identifiability issues
in parameters like scales sigma.

- Bob

[Stefano]

Thanks Bob,

I will articulate better my case with code and images and attached data 

==================================

predicting x (just point)

==================================

data{

  int P;

  int S;

  int R;

  simplex[P] a[S];

  vector[S] y;

  real MU;

  real SIGMA;

}

parameters{

  vector<lower=0.00001>[P] ex;

  real<lower=0.00001> sigma;

  real<lower=0> MU2;

  real<lower=0> SIGMA2;

}

model{

  vector[P] comb;

  

  ex ~ normal(MU, SIGMA);

  sigma ~ cauchy(0, 2.5);

  for(s in 1:S) {

    for(p in 1:P) comb[p] = ex[p] + log(a[s,p]);

    y[s] ~ normal(log_sum_exp(comb), sigma); 

  }

}

The convergence is pretty bad but the results get close to the real

 

PREDICTED  8.417581e+00 6.799527e+00  5.573320e+0  4.938242e+00  5.508769e+00

REAL 8.668656 6.934466   4.307549    4.312484    4.607683

How come the convergence and precision is so low?

==================================

predicting x distribution (meand and STANDARD DEVIATION)

==================================

should i use something like?

 for(s in 1:S) {

  for(p in 1:P){

  comb[p] = log(a[p]) + normal_lcdf(y[s] | ex[p], sigma[p]);

  }

    target += log_sum_exp(comb);

  }

Note: in this case I may have to replace normal_lcdf with lognormal_lcdf and y with exp(y) because the mixture is at lognormal scale. The code number one works because of a logaritmic transformation of the system

This might not be trivial, so thanks for your support

[Bob Carpenter]

You shouldn't have parameters MU2 and SIGMA2 that are
unconstrained and not used.

I'm afraid I don't follow why you are running
log_sum_exp. What is comb[p]?

If it's a mixture model at some level, you can
have problems with identifiability.

For efficiency, you don't want to have to keep applying
log(a[s, p]) --- you can do it once in transformed data
to define a variable log_a where log_a[s, p] = log(a[s, p]).

Not sure why you're then considering normal_lcdf.

What's the model trying to do?

- Bob

> On Nov 10, 2016, at 1:34 AM, Stefano <mangiol...@gmail.com> wrote:
>
> Thanks Bob,
>

> I will articulate better my case with code and images and attached data
>
>
> ==================================
> predicting x (just point)
> ==================================
>
> data{
> int P;
> int S;
> int R;
> simplex[P] a[S];
> vector[S] y;
> real MU;
> real SIGMA;
> }
> parameters{
> vector<lower=0.00001>[P] ex;
> real<lower=0.00001> sigma;
> real<lower=0> MU2;
> real<lower=0> SIGMA2;
> }
> model{
> vector[P] comb;
>
> ex ~ normal(MU, SIGMA);
> sigma ~ cauchy(0, 2.5);
>
> for(s in 1:S) {
> for(p in 1:P) comb[p] = ex[p] + log(a[s,p]);
> y[s] ~ normal(log_sum_exp(comb), sigma);
> }
> }
>
>
>
>
>
>

> <code1.R><data1.RData>

[stefano]

Sorry for the confusion, trying to simplify my model as example was not the best idea.

#1 Let's take stan manual example for linear regression (alpha = 0, beta_1 and beta_2)

--------------------------------

data {
int<lower=0> N;
vector[N] x_1;
vector[N] x_2;
vector[N] y;
}
parameters {
real beta_1;
real beta_2;
real<lower=0> sigma;
}
model {
y ~ normal(beta_1 * x_1 + beta_2 * x_2, sigma);
}
--------------------------------

#2 let's suppose I have many beta as data and one x as parameter. The model becomes

--------------------------------
data {
int<lower=0> N;
vector[N] beta_1;
vector[N] beta_2;
vector[N] y;
}
parameters {
real x_1;
real x_2;
real<lower=0> sigma;
}
model {
y ~ normal(beta_1 * x_1 + beta_2 * x_2, sigma);
}
--------------------------------

#3 Now, let's suppose I don't want just the real x, but I want x_mean and x_sd. How do I do it?

--------------------------------
data {
int<lower=0> N;
vector[N] beta_1;
vector[N] beta_2;
vector[N] y;
}
parameters {
real x_1_mean;
real x_2_mean;
real x_1_sd;
real x_2_sd;
real<lower=0> sigma;
}
model {
y ~ normal(beta_1 * ??? + beta_2 * ???, sigma);
}
--------------------------------

I hope this proof of principle is a bit clearer.

Thanks!

Stefano Mangiola

Data analyst (MSc, MPhil, PhD candidate)
@ Melbourne University

Hovens lab - Papenfuss lab

Royal Melbourne Hospital - WEHI

+61 (0)450135954

> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
> <code1.R><data1.RData>

--
You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.
To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

[Bob Carpenter]

> On Nov 10, 2016, at 7:50 PM, stefano <mangiol...@gmail.com> wrote:
>
> Sorry for the confusion, trying to simplify my model as example was not the best idea.
>
> #1 Let's take stan manual example for linear regression (alpha = 0, beta_1 and beta_2)
>
> --------------------------------
> data {
> int<lower=0> N;
> vector[N] x_1;
> vector[N] x_2;
> vector[N] y;
> }
> parameters {
> real beta_1;
> real beta_2;
> real<lower=0> sigma;
> }
> model {
> y ~ normal(beta_1 * x_1 + beta_2 * x_2, sigma);
> }
> --------------------------------
>
> #2 let's suppose I have many beta as data and one x as parameter.

That doesn't make sense for the linear regression. The
beta are unobserved parameters that are fit by observing
x and y.

> The model becomes
>
> --------------------------------
> data {
> int<lower=0> N;
> vector[N] beta_1;
> vector[N] beta_2;
> vector[N] y;
> }
> parameters {
> real x_1;
> real x_2;
> real<lower=0> sigma;
> }
> model {
> y ~ normal(beta_1 * x_1 + beta_2 * x_2, sigma);
> }
> --------------------------------

That just turns around the names, making beta the data and x_1 the
parameters.

> #3 Now, let's suppose I don't want just the real x, but I want x_mean and x_sd. How do I do it?
>
> --------------------------------
> data {
> int<lower=0> N;
> vector[N] beta_1;
> vector[N] beta_2;
> vector[N] y;
> }
> parameters {
> real x_1_mean;
> real x_2_mean;
> real x_1_sd;
> real x_2_sd;
> real<lower=0> sigma;
> }
> model {
> y ~ normal(beta_1 * ??? + beta_2 * ???, sigma);
> }
> --------------------------------
>
> I hope this proof of principle is a bit clearer.

I'm afraid not.

I'd recommend thinking generatively in terms of how
everything comes about: parameters beta and sigma
from priors, predictors x from priors, and then y from
x, beta, and sigma.

Looked at that way, the density function f(beta, sigma, x, y)
never changes, just what you observe does.

Only the x and y are set up to be directly observable here,
which is why I'm confused about what you're trying to do.

And I don't see where means come into play at all.

- Bob

[stefano]

Thanks for your answer,

let me show you with a graphics I made the generative process

My goal is to estimate 3 distributions from two sources of data a proportion (simplex) and a value (real)

Thanks

Stefano Mangiola

Data analyst (MSc, MPhil, PhD candidate)
@ Melbourne University

Hovens lab - Papenfuss lab

Royal Melbourne Hospital - WEHI

+61 (0)450135954

[Stefano]

Thanks for your answer,

let me show you with a graphics I made the generative process

My goal is to estimate 3 distributions from two sources of data couples proportions (simplex) and values (real)

Thanks

[Bob Carpenter]

Sorry for being dense here, but I'm still not sure I get it.

Is the generative process:

1. draw parameters for components theta[1], theta[2], and theta[3]

2. for data points n in 1:N

a. draw a simplex phi[n] from ???

b. draw from each component
z[n, 1] ~ foo(theta[1])
...
z[n, 3] ~ foo(theta[3])

c. set y[n] = sum(z[n, 1:3] .* phi[n])

In that last case, it's an elementwise product.

And then the inference problem is:

Given observations of phi[1:N] and y[1:N], estimate theta[1:3]?

In other words, it's not a traditional mixture, but rather
we're only observing the sum of three variables, each of which
has a distribution? And did I get the way to produce y right?

- Bob

> On Nov 11, 2016, at 11:15 PM, stefano <mangiol...@gmail.com> wrote:
>
> Thanks for your answer,
>
> let me show you with a graphics I made the generative process
>

> <image.png>

> To unsubscribe from this group and all its topics, 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.
>
>

> --
> 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.

[stefano]

No worries, I'm probably not great in explaining it.

Basically I have 

> To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>
>
> --

> 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+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

--

You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

[stefano]

2016-11-12 16:36 GMT+11:00 stefano <mangiol...@gmail.com>:

No worries, I'm probably not great in explaining it.

2016-11-12 16:30 GMT+11:00 Bob Carpenter <ca...@alias-i.com>:

Sorry for being dense here, but I'm still not sure I get it.

Is the generative process:

1.  draw parameters for components theta[1], theta[2], and theta[3] 

2.  for data points n in 1:N

    a.  draw a simplex phi[n] from ???

    b.  draw from each component
           z[n, 1] ~ foo(theta[1])
           ...
           z[n, 3] ~ foo(theta[3])

    c.  set y[n] = sum(z[n, 1:3] .* phi[n])

In that last case, it's an elementwise product.

And then the inference problem is:

Given observations of phi[1:N] and y[1:N], estimate theta[1:3]?

In other words, it's not a traditional mixture, but rather
we're only observing the sum of three variables, each of which
has a distribution?  And did I get the way to produce y right?

I think you gut it :)

I try to repeat your model in my words, and should match

Basically I have

theta1 drawn from ~ lognormal(mu1, sigma1)

theta2 drawn from ~ lognormal(mu2, sigma2)

theta3 drawn from ~ lognormal(mu3, sigma3)

y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])

The key 6 parameters I want are mu1,2,3 and sigma1,2,3

[Bob Carpenter]

Do you observe y and phi?

- Bob

> > To unsubscribe from this group and all its topics, 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.
> >
> >
> > --
> > 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.
>
> --
> You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

> To unsubscribe from this group and all its topics, 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.
>
>
>
> --
> 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.

[stefano]

Both observed.

[Bob Carpenter]

theta1 drawn from ~ lognormal(mu1, sigma1)
theta2 drawn from ~ lognormal(mu2, sigma2)
theta3 drawn from ~ lognormal(mu3, sigma3)

y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])

If y and phi are both known, then it gets tricky. You
only have two degrees of freedom, and given that the theta
are positive (they're lognormal), you have

0 < theta[n, 1] * phi[n, 1] < y[n]

and then

0 < theta[n, 2] * phi[n, 2] < y[n] - theta[n, 1] * phi[n, 1]

and then

theta3[n] = y[n] - (theta1[n] * phi[n, 1]) + (theta2[n] * phi[n, 2]);

You have to do those two constraints custom, define theta3 as
a transformed parameter, and then do the Jacobian adjustments.

Possible, but a pain to code.

- Bob


> On Nov 12, 2016, at 5:07 AM, stefano <mangiol...@gmail.com> wrote:
>
> Both observed.

[stefano]

2016-11-13 13:26 GMT+11:00 Bob Carpenter <ca...@alias-i.com>:

  theta1 drawn from ~ lognormal(mu1, sigma1)
  theta2 drawn from ~ lognormal(mu2, sigma2)
  theta3 drawn from ~ lognormal(mu3, sigma3)

  y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])

If y and phi are both known, then it gets tricky.  You
only have two degrees of freedom, and given that the theta
are positive (they're lognormal), you have

  0  <  theta[n, 1] * phi[n, 1]  <  y[n]

and then

  0  < theta[n, 2] * phi[n, 2] < y[n] - theta[n, 1] * phi[n, 1]

and then

 theta3[n] = y[n] - (theta1[n] * phi[n, 1]) + (theta2[n] * phi[n, 2]);

You have to do those two constraints custom, define theta3 as
a transformed parameter, and then do the Jacobian adjustments. 

Possible, but a pain to code.

I see your point. 

Let's start with saying that this is a last step of a method that already outperformed a one-year-old Stanford-made alternative (Nature methods), and I will do whatever it takes to complete and publish this high impact advance and paper. 

So if it will be painful to code, well I will do it (maybe with some help :) )

 

Said that I will raise 3 points:

= 1 ==================================================

What if I say that actually, phi is a simplex for which 

y[n] = ( (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3]) ) + theta_background * phi_background

Does this make any difference?

= 2 ==================================================

What if I take a less deterministic approach

instead of  

y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])

I use

y[n] ~ distribution( (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])  , error )

= 3 ==================================================

The number of dimension of phi and theta is not established, could be 3, 4, 5, 6...

=====================================================

Thanks a lot for the support Bob!

Bw

The other thing that  me

- Bob


> On Nov 12, 2016, at 5:07 AM, stefano <mangiol...@gmail.com> wrote:
>
> Both observed.

>
>
> --
> 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+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.
To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

[stefano]

2016-11-13 19:22 GMT+11:00 stefano <mangiol...@gmail.com>:

2016-11-13 13:26 GMT+11:00 Bob Carpenter <ca...@alias-i.com>:

  theta1 drawn from ~ lognormal(mu1, sigma1)
  theta2 drawn from ~ lognormal(mu2, sigma2)
  theta3 drawn from ~ lognormal(mu3, sigma3)

  y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])

If y and phi are both known, then it gets tricky.  You
only have two degrees of freedom, and given that the theta
are positive (they're lognormal), you have

  0  <  theta[n, 1] * phi[n, 1]  <  y[n]

and then

  0  < theta[n, 2] * phi[n, 2] < y[n] - theta[n, 1] * phi[n, 1]

and then

 theta3[n] = y[n] - (theta1[n] * phi[n, 1]) + (theta2[n] * phi[n, 2]);

You have to do those two constraints custom, define theta3 as
a transformed parameter, and then do the Jacobian adjustments. 

Possible, but a pain to code.

I see your point. 

Let's start with saying that this is a last step of a method that already outperformed a one-year-old Stanford-made alternative (Nature methods), and I will do whatever it takes to complete and publish this high impact advance and paper. 

So if it will be painful to code, well I will do it (maybe with some help :) )

 

Said that I will raise 3 points:

= 1 ==================================================

What if I say that actually, phi is a simplex for which 

y[n] = ( (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3]) ) + theta_background * phi_background

sorry I correct 

y[n] = ( (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3]) ) * (1-phi_background) + theta_background * phi_background

[Bob Carpenter]

> On Nov 13, 2016, at 3:36 AM, stefano <mangiol...@gmail.com> wrote:
>
>
> What if I take a less deterministic approach
>
> instead of
>
> y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])
>
> I use
>
> y[n] ~ distribution( (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3]) , error )

That can work especially if you put a small prior on
the error scale or just set it.

- Bob

[Bob Carpenter]

> On Nov 13, 2016, at 3:36 AM, stefano <mangiol...@gmail.com> wrote:
>

That'd just be another way of adding noise, right? And the
phi_background or theta_background would vary with n?

- Bob

[stefano]

Thanks,

It works, now my problem is that I don't want just theta1, but

- theta1_mean

- theta1_sd

 

Any idea on how could I code it? This is the initial bit where I was getting lost

[Bob Carpenter]

Why do you want mean(theta) and sd(theta)? If they come from
a hierarchical normal prior, why not use its mean and standard deviation?

If you really just want those quantities, Stan has functions mean(theta)
and sd(theta) for sample means and standard deviations. So you can define
those as transformed parameters or generated quantities. But I wouldn't
try to put a distribution on them as it will no longer be a clean generative
model. And if you do insist on putting a distribution on them, then you
need to take care of the Jacobian adjustment manually by incrementing the
target log density.

- Bob

> On Nov 13, 2016, at 11:50 PM, stefano <mangiol...@gmail.com> wrote:
>
> 2016-11-14 14:41 GMT+11:00 Bob Carpenter <ca...@alias-i.com>:
>
> > On Nov 13, 2016, at 3:36 AM, stefano <mangiol...@gmail.com> wrote:
> >
> >
> > What if I take a less deterministic approach
> >
> > instead of
> >
> > y[n] = (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3])
> >
> > I use
> >
> > y[n] ~ distribution( (theta1 * phi[n,1]) + (theta2 * phi[n,2]) + (theta3 * phi[n,3]) , error )
>
> That can work especially if you put a small prior on
> the error scale or just set it.
>
> - Bob
>
>
> Thanks,
>
> It works, now my problem is that I don't want just theta1, but
>
> - theta1_mean
> - theta1_sd
>
> Any idea on how could I code it? This is the initial bit where I was getting lost
>
> --
> You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

> To unsubscribe from this group and all its topics, 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.
>
>

> --
> 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.

[stefano]

2016-11-16 5:45 GMT+11:00 Bob Carpenter <ca...@alias-i.com>:

Why do you want mean(theta) and sd(theta)?  If they come from
a hierarchical normal prior, why not use its mean and standard deviation?

You are right, sorry

something like would work?

theta[p] ~ lognormal(mu[p], sigma[p])
...
y[n] ~ distribution( (theta[1] * phi[n,1]) + (theta[2] * phi[n,2]) + (theta[3] * phi[n,3])  , error )

Thanks

 

> To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>
>
> --

> 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+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

--

You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

[Bob Carpenter]

Sure, that'll work. It might have identifiability issues,
but it's well formed both conceptually and for Stan.

- Bob

> > To unsubscribe from this group and all its topics, 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.
> >
> >
> > --
> > 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.
>
> --
> You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

> To unsubscribe from this group and all its topics, 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.
>
>
> --
> 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.

[stefano]

Thanks Bob.

I'm exploring more this route and apparent  problems with vectorization.

Bw.

> > To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

> > To post to this group, send email to stan-...@googlegroups.com.
> > For more options, visit https://groups.google.com/d/optout.
> >
> >
> > --
> > 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+unsubscribe@googlegroups.com.

> > To post to this group, send email to stan-...@googlegroups.com.
> > For more options, visit https://groups.google.com/d/optout.
>
> --
> You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

> To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>
>
> --
> 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+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/em6RroJpUvM/unsubscribe.

To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.