2
$\begingroup$

I have a sequence of binary samples and I want to model how this sequence is changing. The change should reflect whether it's a positive or negative change.

For example, sequence 1:

1 1 1 1 0 1 1 1 0 1 1 1 0 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1

is changing positively since we're confident in appearance of 1's before and become even more confident after some point. However, sequence 2:

1 1 1 1 0 1 1 1 0 1 1 1 0 1 1 0 1 0 1 0 0 0 1 1 0 0 0 0 0 0 0 0

is changing negatively since the number of 0's is increasing, as opposed to our prior confidence in 1. Is there a metric that models this behavior?

I've been looking into sample entropy but it seems to model change of frequency but does not distinguish positive and negative changes. The reason why I'd like to model this difference is I want to know whether there is a positive or negative effect on the underlying graphical model that generates these samples so that I can detect these underlying changes. Thanks!

$\endgroup$
  • $\begingroup$ Are you analyzing a sequence whose length is changing or are you just considering the last, say, $k$ elements? This seems very much like a coin toss problem, where your coin has a stochastic probability $p(t)$. In the easiest case, I think you can just update an estimate of this probability and see where it changes significantly (say, via the Bayes theorem). $\endgroup$ – Néstor Jul 17 '12 at 2:14
  • $\begingroup$ I'm considering the last $k$ elements and want to know whether $p(x)$ is changing based on these samples. I have many random variables and want to map them to a real number so that I can know which random variables are changing most rapidly so that I can generate more samples for them. (these random variables are nodes of a graphical model) $\endgroup$ – Yang Jul 17 '12 at 3:23
  • 1
    $\begingroup$ Possibly related: this question. $\endgroup$ – MånsT Jul 17 '12 at 7:47
1
$\begingroup$

What you are describing is a changepoint detection problem. Your data $y_1,\dots,y_n$ follows a Bernoulli distribution with unknown parameter $p$. You are interested in finding such point $1 \le \tau \le n$ that marks change in the distribution, so that values before it occur with probability $p_1$ and after it, with probability $p_2$. To find it, you can use likelihood ratio test based on log-likelihoods $\ell$,

$$ \underset{\tau}{\mathrm{arg\ max}} \{ \ell(y_{1:\tau-1}) + \ell(y_{\tau:n}) - \ell(y_{1:n}) \} \tag{1}$$

where $\ell$ is a logarithm of a standard binomial likelihood

$$ p^k (1-p)^{n-k} \tag{2}$$

where $p$ is a probability of success (estimated from data), $k$ is number of successes and $n$ is a total number of trials (length of sequence).

Moreover, formula (1) is a likelihood ratio and can be used for hypothesis testing, where $\text{LR} > \beta$ for some threshold $\beta$ can be used for testing if changepoint occurred.

The same in R code:

binchp <- function(x) {
  n <- length(x)
  llik <- function(x) {
    n <- length(x)
    k <- sum(x)
    p <- k/n
    log(p)*k + log(1-p)*(n-k) 
  }
  ll <- vapply(2:(length(x)-1), function(i) llik(x[1:i]) + llik(x[i:n]), numeric(1))
  which.max(ll)
}

Notice that in the code above likelihood of null model (no changepoint) was omitted, since it is the same in all cases. If you conduct simple simulation you'll see that it works very well for such data and most estimated values are close to the true changepoint.

p1 <- 0.5
p2 <- 0.3

x1 <- rbinom(150, 1, p1)
x2 <- rbinom(190, 1, p2)
xx <- c(x1, x2)

res <- replicate(5000, {
  x1 <- rbinom(150, 1, p1)
  x2 <- rbinom(190, 1, p2)
  xx <- c(x1, x2)
  binchp(xx)
})

hist(res-150, 100)

Distribution of errors

This can be easily extended to more then one changepoint, but you would rather need more efficient optimization algorithm then simple brute-force search.

To read more about changepoint models you can check paper by Rebecca Killick and Idris Eckley (2013) describing their changepoint package.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.