Я пытаюсь решить систему из n уравнений в R. n-1 уравнений нелинейны, а последнее - линейно. Если это помогает, это ограниченная задача оптимизации, где уравнение n-1 является условиями первого порядка, а последнее - бюджетным ограничением. Н-1 нелинейные уравнения могут быть охарактеризованы следующим образом: нелинейные уравнения
Если изображение не отображается, его можно определить поэлементно, например:
v_i * epsilon_i * cos (2 / pi * e_i / delta_i) -lambda = 0
, где epsilon, v, e и delta являются векторами n- 1 измерение, и лямбда является скаляром, общим для всех уравнений).
Последнее уравнение является простым линейным уравнением вида | e | = c. То есть сумма всех элементов в e - это некоторая известная c, называемая ниже parms [1,4] или "бюджет".
Я заинтересован в решении для вектора e и постоянной лямбды, относиться ко всему остальному как данному.
Я пытался использовать rootSolve's multi root. Чтобы сделать это, я определяю один вектор X, который должен быть вектором e, с добавлением лямбды к нему, так что multi root решает для x, и получает n уравнений в виде списка. Все параметры сохраняются в матрице под названием «parms»
Сначала я определяю n-1 нелинейных уравнений
convex_focs <- function(x = numeric(),parms = numeric()){
deltas = parms[,1]
v = parms[,2]
lambda = x[1]
e = x[2:length(x)]
epsilon_2 = exp(parms[,3]) - parms[,1]
return(epsilon_2*cos((pi/2)*(e/deltas))-lambda)
}
В этом уравнении используется матричная запись, и она отлично работает сама по себе , Затем я определяю последнее линейное уравнение:
convex_budget <- function(x = numeric(),parms = numeric()){
e = x[2:length(x)]
return(sum(e)-parms[1,4])
}
Затем я пробовал convex_system <- function(x,parms) c(convex_focs,convex_budget )
и вызываю:
multiroot(f = convex_system, maxiter = 500000, start = c(0,rep(budget/length(parms[,1]),length(parms[,1]))), parms = parms[,])
Это, конечно, не работает, так как rootSolve распознает convex_system
как два уравнения, но X как n-мерность.
Если я опущу последнее уравнение и трактую лямбду как заданную (таким образом, только решая нелинейные уравнения) , я могу получить решение . Но, конечно, это нехорошо, потому что я не знаю лямбду.
Итак, мой первый вопрос: 1. Как мне создать список функций из моих векторов, которые rootSolve распознает как систему? Я пытался использовать lapply или виньетки, чтобы создать список convex_focs
уравнений, по одному на каждый элемент в векторе, но не мог придумать, как заставить это работать. 2. Почему он распознал мою исходную функцию convex_focs
как систему уравнений, но когда я добавил convex_budget
, она перестала работать?
Я тогда (в отчаянии ...) попытался вручную определить набор функции, глядя только на 3 нелинейных, а не n-1. Это связано с тем, что мой список функций будет выглядеть как руководство и другие решения, которые я нашел в Интернете :
convex_system <- function(x,parms) c(F1 =
function(x =x,parms = parms){
deltas = parms[1,1]
v = parms[1,2]
lambda = x[1]
e = x[2]
epsilon_2 = exp(parms[1,3]) - parms[1,1]
return(v*epsilon_2*cos((pi/2)*(e/deltas))-lambda)
}
,
F2 =
function(x = x,parms = parms){
deltas = parms[2,1]
v = parms[2,2]
lambda = x[1]
e = x[3]
epsilon_2 = exp(parms[2,3]) - parms[2,1]
return(v*epsilon_2*cos((pi/2)*(e/deltas))-lambda)
}
,
F3 =
function(x = x,parms = parms){
deltas = parms[3,1]
v = parms[3,2]
lambda = x[1]
e = x[4]
epsilon_2 = exp(parms[3,3]) - parms[3,1]
return(v*epsilon_2*cos((pi/2)*(e/deltas))-lambda)
}
,
F_budget = function(x = x,parms = parms){
e = x[2:length(x)]
return(sum(e)-parms[1,4])}
)
И вызов multiroot(f = convex_system, maxiter = 500000, start = c(0,rep(budget/length(parms[1:3,1]),length(parms[1:3,1]))), parms = parms[1:3,])
Когда я запускаю это Я получаю сообщение об ошибке
Ошибка в краже (y, times, fun c, parms = parms, ...): REAL () может быть применена только к 'цифре c' а не «список»
Что я действительно не понимаю - как может список функций не относиться к классу «список»? Итак, мой второй вопрос:
Как создать список функций, если они не являются простыми однострочными функциями (как в приведенных выше ссылках)
Наконец, я очень признателен за любые указания о том, как Лучше решить такие проблемы в R.
Спасибо за любую помощь!