# Time Series sampled at varying frequency - Employ Linear Mixed Model to compare trends?

I hope that this question has not be asked like this elsewhere, if so I could not find it during my google research..

I have the following problem: I have data sampled from different sensors("ID") placed in different machines("machine") across the population. The data is coming in very sparse, so I do not have individual time series for every sensor. The value the sensors measure is continuos ("value" in the example data)

My goal is now, with the data accessible atm, to cluster these sensors (in different product lines) and obtain a time series for every type (e.g. types "A","B","C" below). The objective would then be to analyse whether they expose a different trend by sensor type over time.

The data could look like this:

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd

dates = pd.date_range(start='01-01-2015',periods=2000,freq='D')

IDs = np.arange(1,177)
machines = np.arange(1,21)
np.random.seed(2018)

## just creating some random data that resembles the structure of my data
rand_dates = np.random.choice(dates,200)

## every machine has sensors of all 4 types, a sensor can move from one machine to another
rand_machines = np.array([[np.random.choice(machines)]*4 for i in range(50)])
rand_machines = rand_machines.reshape(200,)
rand_IDs = [np.random.choice(IDs[(IDs>=((i-1)%5)*46)&(IDs<=(i%5)*46)])  for j in range(50) for i in range(1,5)]

rand_values = np.random.rand(200)
df = pd.DataFrame(data=[rand_machines,np.array(rand_IDs),rand_values]).transpose()
df.index = rand_dates
df.columns = ['ID','value','machine']

df.columns = ['machine','ID','value']
df.index.name = 'date'
df['sensor_type'] ='A'
df['sensor_type'][df.ID>46]='B'
df['sensor_type'][df.ID>92]='C'
df['sensor_type'][df.ID>138]='D'

machine     ID     value sensor_type
date
2018-06-28     11.0    7.0  0.209688           A
2016-09-06     11.0   52.0  0.290023           B
2015-08-15     11.0   93.0  0.957608           C
2017-02-16     11.0  154.0  0.576920           D
2017-07-06      1.0   40.0  0.601943           A


One problem (as illustrated by the example-data above) is that the samples are not sampled at any regular frequency, they just come in when they are reported (externally).

I hoped to treat them as time series (i.e. 4 different time series, one for each sensor type) and then try something like DTW with subsamples growing in length. My concern with this was the different sample frequency and the fact that they are basically sampled at different machines, which raised the question whether I can treat this as a time series at all.

Then I came across Linear Mixed Effects (I'm not really familiar with it, so I do not know the whole statistics behind and I do not know if I will have the time to read myself completely into it, as the problem described here is not of the highest importance).

Anyway, I was surprised by the possibility to include almost anything into the linear mixed models, as people seem to be doing, so I gave it a try in Python. Remember, my goal is to identify whether the subtypes expose different trends over time - so i thought I could include the time as an independent variable by arbitralily converting it to seconds:

df = df.reset_index(drop=False)
df.date = (df.date - min(df.date))/np.timedelta64(1, 's')


Then I would include this as an indepent variable in the Linear Effects Model and check for the interaction with sensor_type (grouped by machine, as every machine contains all different sensor types). This is how I would specify the model:

import statsmodels.api as sm
import statsmodels.formula.api as smf

md = smf.mixedlm("value ~ date*sensor_type",data=df,groups = df["machine"])
mdf = md.fit()
print(mdf.summary())

Mixed Linear Model Regression Results
================================================================
Model:               MixedLM    Dependent Variable:    value
No. Observations:    200        Method:                REML
No. Groups:          18         Scale:                 0.0798
Min. group size:     4          Likelihood:            -116.5445
Max. group size:     24         Converged:             No
Mean group size:     11.1
----------------------------------------------------------------
Coef.  Std.Err.   z    P>|z| [0.025 0.975]
----------------------------------------------------------------
Intercept              0.525    0.071  7.441 0.000  0.387  0.664
sensor_type[T.B]       0.087    0.114  0.765 0.444 -0.136  0.309
sensor_type[T.C]      -0.016    0.109 -0.150 0.881 -0.230  0.198
sensor_type[T.D]      -0.048    0.111 -0.437 0.662 -0.265  0.168
date                  -0.000    0.000 -0.104 0.917 -0.000  0.000
date:sensor_type[T.B] -0.000    0.000 -0.378 0.705 -0.000  0.000
date:sensor_type[T.C]  0.000    0.000  0.987 0.323 -0.000  0.000
date:sensor_type[T.D]  0.000    0.000  0.171 0.864 -0.000  0.000
groups RE              0.000    0.008
================================================================


This does not converge, but I guess it is because of the randomness of the data (on my data it does converge, the structure of the data is similiar)

Now I Have some questions:

1) Does this make sense what I am doing above, i.e. am I adressing the question I want to adress (different trends over time for different sensor types?)

2) How could I now include additional random effects in this model? I think the sampling-time should also be somehow modelled as having a random effect ( as I have no possibility to control it), but how could I do this?

Also I have the sensor IDs as another group that is nested in the machine, but at various sample points this can change (sometimes sensor gets used in another machine). So can I include this via the "vcf" parameter? (when I do my model does not converge...

3) Is there any requirements to the number of samples per grouping variable (e.g. do I need at least X samples from a machine to make it a meaningful group)

4) Reporting the effects: Can I use the typical p-value heuristics here to report a significant effect? Or how should I determine whether there is a different trend or not?

Thanks a lot in advance to you guys ! (I can use preferably Python or also R for my analysis)

I would be also glad for any comprehensive links how to deal with similiar kind of problems. Thank you guys!

• If the data exhibits strong dependance on time, perhaps you could also consider fitting a model on all of the data, disregarding the sensor type, that looks vaguely like: $$value = long\_term\_trend + seasonal\_variation + error$$ Then, after fitting a TS model on all the data, you could refit the same model on subsets of the data corresponding to each individual sensor type, and see if the parameter estimates differ significantly from the assumption that they all produce similar estimates for value. A linear model could be reasonable, if the time factor isn't a huge deal or is linear. – InfProbSciX Oct 12 '18 at 14:13
• Thanks for the idea! But this still leaves the problem of unequal sampling rate, doesn't it? The approaches I found for dealing with unequal frequency e.g. for an ARIMA depend on interpolating to have an equal frequency... but I really have gaps varying from hours to multiple days, so I don't know whether I could interpolate without the majority of my data then being interpolations. Do you know any models that can do the decomposition you specified, without relying on regular sample rate? (besides, I also have 4 measures per machine on a single time point) – a_student Oct 19 '18 at 10:12
• I do not assume a regular sampling rate - if you can assume that different machines of the same sensor type sampled at the same time produce similar results (which could be verified by looking at the data points where the gaps were really small), then you could simply use a periodogram (spectrum in R) to identify the dominant frequencies to identify seasonality. The trend and error terms don't depend on the unequal sampling anyway. If there is a significant amount of error between the machines, it'd be harder to identify the frequencies, you'd have to smooth the series first somehow. – InfProbSciX Oct 23 '18 at 9:23
• Another note, the kind of seasonal variations I'm suggesting are perhaps parametric terms like $sin(2*pi*f*x + c)$ or a non-parametric estimator like a GP with an exponential-sine kernel. For the errors, you could still use an ARIMA model if you so desired, by treating the missing values as missing (in PyMC, BUGS, JAGS or Stan for example) – InfProbSciX Oct 23 '18 at 9:28