How to infer size of a rolling window using Stan?

Hi everyone, I’m wondering if it’s possible to use Stan for the following problem:
I have 2 time-series y(t) and z(t) and the relation I want to model is

z(t) = \int_{t - b}^{t} y(s)\mathrm{d}s ,

with b a positive real number. The idea I’m using now is basically re-writing the integral as a dot product of 2 vectors: z(t) = v(t, b)\cdot y(t) where the vector v(t, b) equals to 1 in the interval (t - b, t) and then it decays smoothly to zero outside of those values. Does that sound like a reasonable approach? I feel like there must be a better way so I’m looking for suggestions (below my current code).

functions {
    real smoothstep(real x) {
        // A smooth version of the step function
        if (x < 0) {return 0;}
        if (x > 1) {return 1;}
        return 6*x^5 - 15*x^4 + 10*x^3;
    }
    
    vector smoothstep_vector(int T, real x) {
        vector[T] result;
        for (t in 1:T) {
            result[t] = smoothstep(t - x);
        }
        return result;
    }
}

data {
    int T;
    int S;   // ignore the first S points
    int y[T];
    int z[T];
}

transformed data {
    vector[T] vy = to_vector(y);
}

parameters {
    real<lower=0> win;
}

model {
    win ~ normal(10, 1);
    for (t in S:T) {
        z[t] ~ poisson(1e-5 + vy[:t]' * smoothstep_vector(t, t - win));
    }
}

Thank you all! :)

1 Like

This sounds like a relatively difficult problem.

You can probably write this model as a discrete cut point. Check here for an example: 7.2 Change Point Models | Stan User’s Guide

If you don’t have any continuous parameters to sample though, you might not need Stan.

A nitpicky thing, but you’ve written the relationship between z and y as an integral in Math but as a probabilistic thing in the model. Think of this relationship in both cases as probabilistic. Don’t wanna overlook the importance of the likelihood there.

1 Like

is this some kind of (blind)deconvolution problem?

Thank you for the tip, I’ll keep that in mind :)

It’s not immediately clear to me how to write it as a discrete cut point but I’ll think about it, I didn’t know I could use that ternary operator in the model block (it feel somewhat non-differentiable)

It’s possible with if statements and for loops to accidentally write non-differentiable stuff. It can be difficult to debug so you just gotta watch out :P (and void branching on things that are functions of parameters)