I would like to find max value (response) given the size of difference (independent variables).
Here is some data:
x <- "A B C 1 0.63 0.67 0.61 2 0.62 0.64 0.60 3 0.64 0.65 0.59 4 0.70 0.70 0.63 5 0.71 0.73 0.68 6 0.70 0.75 0.69 7 0.71 0.75 0.70 8 0.74 0.76 0.71 9 0.79 0.81 0.74 10 0.80 0.83 0.76 11 0.82 0.84 0.78 12 0.82 0.84 0.80 13 0.83 0.85 0.81 14 0.81 0.88 0.80 15 0.78 0.84 0.77 16 0.75 0.79 0.74 17 0.73 0.77 0.72 18 0.72 0.75 0.71 19 0.73 0.75 0.71 20 0.73 0.75 0.71 21 0.74 0.76 0.72 22 0.72 0.76 0.71 23 0.71 0.74 0.69 24 0.73 0.75 0.70 25 0.78 0.79 0.71 26 0.82 0.84 0.77 27 0.80 0.84 0.78 28 0.77 0.81 0.76 29 0.79 0.81 0.75 30 0.83 0.84 0.78 31 0.86 0.87 0.82 32 0.85 0.87 0.83 33 0.83 0.84 0.82 34 0.78 0.85 0.77 35 0.74 0.80 0.72 36 0.72 0.76 0.71 37 0.74 0.77 0.70 38 0.75 0.75 0.70 39 0.78 0.81 0.72 40 0.78 0.82 0.75" # Or generate it like this x <- data.frame( A = c(0.63, 0.62, 0.64, 0.7, 0.71, 0.7, 0.71, 0.74, 0.79, 0.8, 0.82, 0.82, 0.83, 0.81, 0.78, 0.75, 0.73, 0.72, 0.73, 0.73, 0.74, 0.72, 0.71, 0.73, 0.78, 0.82, 0.8, 0.77, 0.79, 0.83, 0.86, 0.85, 0.83, 0.78, 0.74, 0.72, 0.74, 0.75, 0.78, 0.78), B = c(0.67, 0.64, 0.65, 0.7, 0.73, 0.75, 0.75, 0.76, 0.81, 0.83, 0.84, 0.84, 0.85, 0.88, 0.84, 0.79, 0.77, 0.75, 0.75, 0.75, 0.76, 0.76, 0.74, 0.75, 0.79, 0.84, 0.84, 0.81, 0.81, 0.84, 0.87, 0.87, 0.84, 0.85, 0.8, 0.76, 0.77, 0.75, 0.81, 0.82), C = c(0.61, 0.6, 0.59, 0.63, 0.68, 0.69, 0.7, 0.71, 0.74, 0.76, 0.78, 0.8, 0.81, 0.8, 0.77, 0.74, 0.72, 0.71, 0.71, 0.71, 0.72, 0.71, 0.69, 0.7, 0.71, 0.77, 0.78, 0.76, 0.75, 0.78, 0.82, 0.83, 0.82, 0.77, 0.72, 0.71, 0.7, 0.7, 0.72, 0.75)) And here some adjustment:
data <- read.table(text=x, header = TRUE) data$diff_AC <- with(data, (A-C)) data$diff_AB <- with(data, (A-B)) with(data, plot(A~1, col=1)) with(data, points(B~1, col=2)) with(data, points(C~1, col=3)) Calculate return:
data$retA <- with(data, as.numeric(c(0,diff(A))/lag(A,1))) Now using optim find the size of difference of A vs B and A vs C at which the return of A (retA) is highest given all data of A.
This should be done separately for negative return and positive return.
I have tried this, but I'm not sure how to bring in the return A part into optim.
max.rss <- function(data, par) { with(data, -sum((par[1] * (B - A) + (C - A))^2)) } result <- optim(par = 0, max.rss, data = data, method = "Brent", lower = 0, upper = 1) EDIT:
So the question is at which level (the size of difference) diff_AB and diff_AC should be so that the retA (return of series A) is highest (maximised) and what is the difference size diff_AB and diff_AC at which the retA is at its minimum (the highest negative rate).
with(data, plot(retA ~ diff_AB, ylim=c(-0.1,0.1), xlim=c(-.1,.1))) with(data, points(retA ~ diff_AC, col="red3")) EDIT2:
It's very likely that the question is not well defined, or with current explanation it does not necessarily makes sense.
Any proposition towards analysis/modeling the underlying dynamics is welcome!
EDIT 3:
Here is possible solution somewhat based on existing answer:
data$rank_min <- with(data, ave(retA, diff_AB, FUN=function(x) rank(x, ties.method="min"))) data$rank_max <- with(data, ave(retA, diff_AB, FUN=function(x) rank(x, ties.method="max"))) with(data, data[rank_min==min(rank_min), ]) with(data, data[rank_max==max(rank_max), ]) or (but I'm not sure whether this is fully correct)
diff_binAB <- with(data, unique(diff_AB)) mse <- numeric(length(diff_binAB)) for(i in 1:length(diff_binAB)){ pwise <- with(data, lm(retA ~ diff_AB*(diff_AB < diff_binAB[i]) + diff_AB*(diff_AB >= diff_binAB[i]))) mse[i] <- summary(pwise)[6] } mse <- as.numeric(mse) mse diff_binAB[which(mse==min(mse))] # -0.07 1 Answers
Answers 1
If you don't want to rely on a single observation (like with data[which.max(retA), ]) you can either use 10% quantiles :
with(data, summary(diff_AB[retA < quantile(data$retA, 0.1)])) with(data, summary(diff_AC[retA < quantile(data$retA, 0.1)])) with(data, summary(diff_AB[retA > quantile(data$retA, 0.9)])) with(data, summary(diff_AC[retA > quantile(data$retA, 0.9)])) or estimate a model of the dependence
m1 <- with(data, lm(retA~diff_AC+diff_AB)) summary(m1) m2 <- with(data, lm(retA ~ I(diff_AC+diff_AB))) summary(m2) But you might want to account for the autocorrelation of your time-series residuals:
plot(m2$residuals, type = "l") 
For example, let's assume it is an AR(1) process:
library(astsa) m3 <- with(data, sarima(retA, 1,0,0, xreg=cbind(diff_AB, diff_AC))) 
and then predict retA based on this model.
0 comments:
Post a Comment