Two equations model with error on data

Modeling
1

Hello,

I am fitting a model with two linear equations with shared parameters:

y1 = P1 + P2f(x1) + P3f(x1,x2) + P4f(x2)
y2 = P3 + P4f(x1) + P5f(x1,x2) + P6f(x2)

where (y1, y2) are data with error and (x1, x2) are date with no error. I want to determine parameters P.

What is the best method to include the errors in the data and to match the two equations?

Thank you very much

2

Since you have f(x1,x2), I’ll assume that x1/x2/y1/y2 all have the same number of elements. A not-very-efficient-but-hopefully-pedagogically-illuminating model would be:

data{
  int N ;
  vector[N] y1 ;
  vector[N] y2 ;
  vector[N] fx1 ;
  vector[N] fx2 ;
  vector[N] fx1x2 ;
}
parameters{
  vector[6] p ;
  vector<lower=0>[2] y_noise ;
}
model{
  p ~ ... ; //prior on parameters goes here
  y_noise ~ ... ; //prior on noise magnitude goes here
  y1 ~ normal(
    p[1] + p[2]*fx1 + p[3]*fx1x2 + p[4]*fx2
    , y_noise[1]
  ) ;
  y2 ~ normal(
    p[3] + p[4]*fx1 + p[5]*fx1x2 + p[6]*fx2
    , y_noise[2]
  )
}
1 Like
4

Thank you very much for your answer.

The problem is that I don’t have the same number of elements. I’ve got ~500 groups of (y1/x1/x2) and other ~500 groups of (y2/x1/x2). Something like:

y1 y2 x1 x2
50 NA 30 60
55 NA 20 70
NA 60 30 60
NA 63 30 60

I know how to determine parameter P by each function (y1 and y2) separately, but I get different values of common Ps. The dispersion is big.

5

OK, then you mean you have separate x1 & x2 vectors for each of y1 and y2? In that case you would do (changing some nomenclature for clarity:)

data{
  int N_a ;  
  vector[N_a] y_a ;
  vector[N_a] fx1_a ;
  vector[N_a] fx2_a ;
  vector[N_a] fx1x2_a ;
  int N_b ;
  vector[N_b] y_b ;
  vector[N_b] fx1_b ;
  vector[N_b] fx2_b ;
  vector[N_b] fx1x2_b ;
}
parameters{
  vector[6] p ;
  vector<lower=0>[2] y_noise ;
}
model{
  p ~ ... ; //prior on parameters goes here
  y_noise ~ ... ; //prior on noise magnitude goes here
  y_a ~ normal(
    p[1] + p[2]*fx1_a + p[3]*fx1x2_a + p[4]*fx2_a
    , y_noise[1]
  ) ;
  y_b ~ normal(
    p[3] + p[4]*fx1_b + p[5]*fx1x2_b + p[6]*fx2_b
    , y_noise[2]
  )
}
2 Likes
6

Yes, you got it. I have separate x1 & x2 vectors for each of y1 and y2, as you say. I think I need some time to understand completely the model (sorry, I’m noob at this stuff), but it sounds nice. I will tell you if it works.

Thank you very much.

7

Thank you very much, it’s taking form. But I still have a doubt. The values of the errors of y (what you called y_noise) are data and I want to implement them for the MCMC sampling. So, they shouldn’t be in the parameters, isn’t it? How can I apply it?

8

Oh, if you know the scale of measurement noise, then you can pass those in as data instead of as parameters, yes.

1 Like
9

I got it. It works perfectly.
Thank you very much!