2016-06-02 2 views
0

Я использовал nleqslv пакет, чтобы решить систему двух нелинейных уравненийnleqslv понять логику итерационной процедуры

library(nleqslv) 
read.table("C:\\Users\\Livia\\Desktop\\firm.txt",header=TRUE) 
firm.txt<-"Q K r X.E t E0 
" 
firm<-read.table(textConnection(firm.txt),header=TRUE,stringsAsFactors=FALSE) 
firm 
#Loop with solver for 9 dates: 
for(i in 1:9){ 
K<-firm[i,2] 
r<-firm[i,3] 
X.E<-firm[i,4] 
t<-firm[i,5] 
E0<-firm[i,6] 
BS<-function(x){ 
f<-rep(NA,length(x)) 
f[1]<-(x[1]*pnorm((log(x[1]/K)+(r+(x[2]^2/2))*t)/(x[2]*sqrt(t)))-exp(-r*t)*K*pnorm((log(x[1]/K)+(r+(x[2]^2/2))*t)/(x[2]*sqrt(t))-x[2]*sqrt(t)))-E0 
f[2]<-X.E*E0-pnorm((log(x[1]/K)+(r+(x[2]^2/2))*t)/(x[2]*sqrt(t)))*x[2]*x[1] 
f 
} 
p0<-c((E0+K),X.E*E0/(E0+K)) 
print(p0) 
ans<-nleqslv(p0,BS,control=list(trace=1,allowSingular=TRUE)) 
print(ans) 
} 

Мне бы очень хотелось, чтобы понять логику позади него. У меня есть два неизвестных, и я установить их начальные значения в соответствии с

p0<-c((E0+K),X.E*E0/(E0+K)) 

Как итерация процедура работает? Я уже пытался получить его из пакетов CRAN «nleqslv» 3.0 .pdf, но безуспешно.

+0

Привет Livia, добро пожаловать в StackOverflow. 2 вещи - 1. Вы можете видеть, как эта функция работает, глядя на исходный код, 2. Если вам нужна дополнительная помощь, это на самом деле вопрос для нашего сайта-брака, Cross Validated. stats.stackexchange.com Исходный код доступен для CRAN, Github и запускается функция в R без '()' –

+0

Используемые алгоритмы не документируются в .pdf. Вы можете найти больше в книге Денниса и Шнабеля. Кстати, вы представили эти уравнения в [r-nleqslv-loop-bad-results] (http://stackoverflow.com/questions/36666061/r-nleqslv-loop-bad-results). Почему вы не использовали временные переменные в моем ответе и почему не масштабировали 'f [1]' с 'K'? – Bhas

+0

Привет, Bhas, я не получал никаких уведомлений, поэтому я никогда не видел вашего ответа! Тем не менее, мой вопрос более теоретический, я просто хочу понять механизм: систему двух уравнений с двумя unkowns; задайте начальные значения этих неизвестных, чтобы начать процедуру итерации ... и затем, где эти значения подставляются, чтобы получить все больше и больше результатов, пока они не сближаются? –

ответ

0

См. Мой комментарий выше для caveates.

Это, как говорится, вот то, что под капотом происходит при выполнении этой функции:

> nleqslv 
function (x, fn, jac = NULL, ..., method = c("Broyden", "Newton"), 
    global = c("dbldog", "pwldog", "cline", "qline", "gline", 
     "hook", "none"), xscalm = c("fixed", "auto"), jacobian = FALSE, 
    control = list()) 
{ 
    fn1 <- function(par) fn(par, ...) 
    if (is.null(jac)) { 
     jacfunc <- NULL 
    } 
    else { 
     if (!is.function(jac)) 
      stop("argument 'jac' is not a function!") 
     jacfunc <- function(par) jac(par, ...) 
    } 
    method <- match.arg(method) 
    global <- match.arg(global) 
    xscalm <- match.arg(xscalm) 
    con <- list(ftol = 1e-08, xtol = 1e-08, btol = 0.001, stepmax = -1, 
     delta = "newton", sigma = 0.5, scalex = rep(1, length(x)), 
     maxit = 150, trace = 0, chkjac = FALSE, cndtol = 1e-12, 
     allowSingular = FALSE, dsub = -1L, dsuper = -1L) 
    if (global == "none") 
     con$maxit = 20 
    if (length(control)) { 
     namc <- names(control) 
     if (!is.list(control) || is.null(namc)) 
      stop("'control' argument must be a named list") 
     if (!all(namc %in% names(con))) 
      stop("unknown names in control: ", paste(sQuote(namc[!(namc %in% 
       names(con))]), collapse = ", ")) 
     con[namc] <- control 
    } 
    tmp <- con[["delta"]] 
    if (is.character(tmp)) { 
     k <- match(tolower(tmp), c("cauchy", "newton")) 
     con[["delta"]] <- as.numeric(-k) 
    } 
    else if (!is.numeric(tmp)) 
     stop("control[\"delta\"] should be either character or numeric") 
    on.exit(.C("deactivatenleq", PACKAGE = "nleqslv")) 
    out <- .Call("nleqslv", x, fn1, jacfunc, method, global, 
     xscalm, jacobian, con, new.env(), PACKAGE = "nleqslv") 
    out 
} 
<environment: namespace:nleqslv> 

Если вы хотите, чтобы посмотреть, как другие функции работают, вы можете исследовать их таким же образом (тип имя функции и нажмите ввод без добавления ()).

Если вы хотите сделать более глубокое погружение, полный исходный код здесь:

https://cran.r-project.org/src/contrib/nleqslv_3.0.1.tar.gz

и здесь:

https://github.com/cran/nleqslv

+0

Привет, Hack-R спасибо за ваш ответ. Однако мой вопрос является более фундаментальным, я просто хочу понять механизм: систему двух уравнений с двумя unkowns; задайте начальные значения этих неизвестных, чтобы начать процедуру итерации ... и затем, где эти значения подставляются, чтобы получить все больше и больше результатов, пока они не сближаются? –

+0

@ LiviaG Я вижу. Это определенно вопрос для stats.stackexchange.com, а не для StackOverflow. StackOverflow - это сайт, когда у вас есть ошибки в коде, и вы пытаетесь выяснить, как исправить ваш код. –

Смежные вопросы