A simple, though somewhat inflexible, way would be to
regress $\bar{P}$ on the $I$ series only (no constant). This will minimise
squared differences instead of absolute ones, though.
R example; I start with creating random data:
nobs <- 250 ## length of series
ns <- 3 ## number of I series
P <- c(1, cumprod(1 + rnorm(nobs, sd = 0.01)))
M <- sum(range(P))/2 ## midpoint of range
I <- array(rnorm((nobs)*ns, sd = 0.01),
dim = c(nobs, ns))
I <- apply(I, 2, function(x) cumprod(1+x))
I <- rbind(1, I)
plot(P, ylim = range(P, I), type = "l",
col = "darkgreen", lwd = 2)
abline(h = M)
for (i in 1:ns)
lines(I[,i], col = grey(0.7))
## regression
res1 <- lm(rep(M, nrow(I)) ~ -1 + I)
lines(I %*% coef(res1),
col = "blue", lwd = 2,
type = "l")
An alternative way would be to use a generic solver.
Here is an example with Differential Evolution, as
implemented in the R package NMOF
, which I maintain.
## Differential Evolution
library("NMOF")
diff_mean <- function(b, M, I)
sum(abs(M - I %*% b))
res2 <- DEopt(diff_mean,
list(min = rep(-1, ns),
max = rep(1, ns)),
M = M, I = I)
lines(I %*% res2$xbest,
col = "red", lwd = 2)
The advantage of such a solver is that it is much more
flexible: you may use another function to measure the
similarity of the series, or add constraints.