DCC模型问题

epoh 发表于 2012-11-25 13:22
or plot(fit.1)

Make a plot selection (or 0 to exit):

我现在又碰到一个问题了，也是mgarch包的问题。这里面做go-garch是用的基于独立成分ICA的，那是不是就相当于ICA-GARCH啊？那我要做基于主成分的PCA-GARCH应该怎么办？谢谢啊~

july0710 发表于 2012-11-26 11:23
我现在又碰到一个问题了，也是mgarch包的问题。这里面做go-garch是用的基于独立成分ICA的，那是不是就相当 ...

你的PC-GARCH可能是指 Patrick Burns (2005)
  Multivariate GARCH with Only Univariate Estimation
   Multivariate GARCH with Only Univariate Estimation.pdf (57.66 KB)

2012-11-26 14:22:33 上传

那可能要照着他的The PC-GARCH Algorithm
七个步骤做了

epoh 发表于 2012-11-25 13:21
install.packages("rmgarch", repos="http://R-Forge.R-project.org")

Newtons<-function (fun, x, ep=1e-5, it_max=100){
    index<-0; k<-1
    while (k<=it_max){
        x1 <- x; obj <- fun(x);
        x  <- x - solve(obj$J, obj$f);
        norm <- sqrt((x-x1) %*% (x-x1))
        if (norm<ep){
            index<-1; break
        }
        k<-k+1
    }
    obj <- fun(x);
    list(root=x, it=k, index=index, FunVal= obj$f)
}

epoh老师，您好！
       x <- x - solve(obj$J, obj$f);
这一行中，solve是什么意思？
obj$J和obj$f是什么意思？J是Jacobi矩阵，是如何求出Jacobi矩阵的？
在 x <- x - solve(obj$J, obj$f);
x - solve(obj$J, obj$f)为什么要赋给x？
x(k+1)<-x(K)-J(逆）*f(x(K))为什么可以x <- x - solve(obj$J, obj$f)这样写呢？

epoh 发表于 2012-11-26 14:22
你的PC-GARCH可能是指 Patrick Burns (2005)
  Multivariate GARCH with Only Univariate Estimation
...

我再看看

zhangtao 发表于 2012-11-26 15:25
Newtons

Newtons<-function (fun, x, ep=1e-5, it_max=100){
    index<-0; k<-1
    while (k<=it_max){
        x1 <- x; obj <- fun(x);
        x  <- x - solve(obj$J, obj$f);
        norm <- sqrt((x-x1) %*% (x-x1))
        if (norm<ep){
            index<-1; break
        }
        k<-k+1
    }
    obj <- fun(x);
    list(root=x, it=k, index=index, FunVal= obj$f)
}

moment <-function(p){
f<-c(p[1]*p[2]-M1,p[1]*p[2]-p[1]*p[2]^2-M2)
J<-matrix(c(p[2],p[1],p[2]-p[2]^2,p[1]-2*p[1]*p[2]),nrow=2,byrow=T)#Jacobi
list(f=f,J=J)
  }

#####Binomial distribution B(size,prob)
set.seed(123)
y<-rbinom(1000, 20, 0.5);
n<-length(y)
M1<-mean(y);   #9.978
M2<-(n-1)/n*var(y)   #4.977516
p<-c(10,0.7)  # p[1]=k , p[2]=p

Newtons(moment, p)
$root
[1] 19.9101695  0.5011509
$it
[1] 5
$index
[1] 1
$FunVal
[1]  1.826095e-12 -4.090950e-12

##########
答复你的问题
1.obj$J和obj$f是什么意思？如何求出Jacobi矩阵的？
  obj <- fun(x);这里的fun(x)就是moment(p)
  function moment()就已经定义了
  f<-c(p[1]*p[2]-M1,p[1]*p[2]-p[1]*p[2]^2-M2)  #FunVal
  J<-matrix(c(p[2],p[1],p[2]-p[2]^2,p[1]-2*p[1]*p[2]),nrow=2,byrow=T)#Jacobi

2.solve是什么意思？
  solve Solve a System of Equations
      solves the equation a %*% x = b for x, where b can be either a vector or a matrix.
  详细请看帮助?solve
  
3.x - solve(obj$J, obj$f)为什么要赋给x？
  迭代法:
  k=1,x(start value)=[10,0.7]
         经Newton's method,
         x <- x - solve(obj$J, obj$f);
         x=[18.3035020 , 0.4165549]
  这就形成k=2,x(start value)=[18.3035020 , 0.4165549]
                k=3,x(start value)=[19.089003 , 0.527265]
                ...
              当k=5,x(start value)=[19.9101655 , 0.5011505]
                 x value 是[19.9101695  0.5011509]
                其与[19.9101655 , 0.5011505]的norm(4.05753e-06)
                已经小于ep=1e-5,故程序中止

epoh 发表于 2012-11-26 19:24
Newtons

epoh老师，您好！
     您说的我明白了，非常感谢！
在以下程序中，书上说是递归调用，也就是fun 调用fun1,
我不明白：1、为什么要递归调用;
2如何实现fun调用fun1的？
epoh老师，希望您能给我讲讲！
非常感谢！
area <- function(f, a, b, eps = 1.0e-06, lim = 10) {
   fun1 <- function(f, a, b, fa, fb, a0, eps, lim, fun) {
      d <- (a + b)/2; h <- (b - a)/4; fd <- f(d)
      a1 <- h * (fa + fd); a2 <- h * (fd + fb)
      if(abs(a0 - a1 - a2) < eps || lim == 0)
         return(a1 + a2)
      else {
         return(fun(f, a, d, fa, fd, a1, eps, lim - 1, fun)
             + fun(f, d, b, fd, fb, a2, eps, lim - 1, fun))
      }
   }
   fa <- f(a); fb <- f(b); a0 <- ((fa + fb) * (b - a))/2
   fun1(f, a, b, fa, fb, a0, eps, lim, fun1)
}

f <- function(x) 1/x
quad<-area(f,1,5); quad

另外，epoh老师，以下错误如何纠正？
>> A=['first';'second']
??? Error using ==> vertcat
All rows in the bracketed expression must have the same
number of columns.

>>

zhangtao 发表于 2012-11-28 15:50
epoh老师，您好！
     您说的我明白了，非常感谢！
在以下程序中，书上说是递归调用，也就是fun 调用f ...

这个不太好解释
也许你看R function area() Adaptive Numerical Integration,会更清楚
Integrate a function of one variable over a finite range using a recursive adaptive method.
  
The method divides the interval in two and compares the values given by Simpson's rule and the trapezium rule.
If these are within eps of each other the Simpson's rule result is given,
otherwise the process is applied separately to each half of the interval and the results added together.
     http://127.0.0.1:24867/library/MASS/html/area.html

能否把'first','second',列出来以方便解释   

epoh 发表于 2012-11-29 19:11
这个不太好解释
也许你看R function area() Adaptive Numerical Integration,会更清楚
Integrate a fun ...

> library(RODBC)
>   z <- odbcConnectExcel('D:/data/mydata.xls')
>  foo<-sqlFetch(z,'Sheet1')
>  close(z)
>   foo
   year gdp ind city  con
1  2001  21  22    22  22
2  2002  22  23    28  24
3  2003  23  24    34  26
4  2004  24  25    40  28
5  2005  25  26    46  30
6  2006  26  27    52  32
7  2007  27  28    58  34
8  2008  28  29    64  36
9  2009  29  30    70  38
10 2010  30  31    76  40
11 2011  31  32    82  42
12 2012  32  33    88  44
13 2013  33  34    94  46
14 2014  34  35   100  48
15 2015  35  36   106  50
16 2016  36  37   112  52
> model1=lm(gdp~ind+city+con)
错误于eval(expr, envir, enclos) : 找不到对象'gdp'
> model1=lm(gdp~ind+city+con,data=mydata)
错误于is.data.frame(data) : 找不到对象'mydata'
> model1=lm(gdp~ind+city+con,data=foo)
错误于eval(expr, envir, enclos) : 找不到对象'city'
> attach(foo)
> model1=lm(gdp~ind+city+con)
错误于eval(expr, envir, enclos) : 找不到对象'city'
> city
错误: 找不到对象'city'
> foo%city
错误: 意外的输入在"foo%city"里
> fix(foo)
> model1=lm(gdp~ind+con,data=foo)
> moedl1
错误: 找不到对象'moedl1'
> model1

Call:
lm(formula = gdp ~ ind + con, data = foo)

Coefficients:
(Intercept)          ind          con  
         -1            1           NA  

>
epoh老师，您好！
以上city 明明读进去了，怎么会提示找不到呢？

zhangtao 发表于 2012-12-1 15:48
> library(RODBC)
>   z   foo  close(z)
>   foo

foofoofoo.xls (13.5 KB)

2012-12-1 20:52:23 上传

library(RODBC)
channel <- odbcConnectExcel("foo.xls")
sqlTables(channel)
#                       TABLE_CAT TABLE_SCHEM TABLE_NAME   TABLE_TYPE REMARKS
# C:\\Users\\...\\Documents\\foo        <NA>       foo$ SYSTEM TABLE    <NA>

foo=sqlFetch(channel, "foo")
foo
#  year        gdp        ind      city         con
#1 2000 0.08431279 0.09429051 0.3622000 0.006607202
#2 2001 0.08300318 0.08441516 0.3766000 0.007220217
#3 2002 0.09082068 0.09829292 0.3908978 0.001303356
#4 2003 0.10025379 0.12671907 0.4053023 0.011064107
#5 2004 0.10085040 0.11112455 0.4176001 0.015448986
#6 2005 0.11310035 0.12082325 0.4299000 0.012044374
#7 2006 0.12676534 0.13391125 0.4434287 0.013780144
#8 2007 0.14162395 0.15063330 0.4588900 0.010812481
#9 2008 0.09634668 0.09877462 0.4698913 0.010299511

model1=lm(gdp~ind+city+con,data=foo)
model1

Call:
lm(formula = gdp ~ ind + city + con, data = foo)

Coefficients:
(Intercept)          ind         city          con  
   -0.03515      0.75016      0.13689     -0.28188