File: grad_fun.R

package info (click to toggle)
r-cran-regsem 1.6.2+dfsg-1
  • links: PTS, VCS
  • area: main
  • in suites: bullseye, sid
  • size: 496 kB
  • sloc: cpp: 263; ansic: 15; makefile: 2
file content (96 lines) | stat: -rw-r--r-- 3,150 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96

grad_fun = function(par,ImpCov,SampCov,Areg,Sreg,A,S,F,
                    A_fixed,A_est,S_fixed,S_est,lambda,type,pars_pen,I){


  # m = dim(ImpCov)[1]
  #grad = rep(0,length(par))
  #facLoad <- matrix(c(1,par[1:(nload-1)]),ncol(A)-1,nfac) #extract just the factor loadings
  #facLoad <- A[1:nload,ncol(A)]
  #facLoad <- A[1:nvar,(ncol(A)-nfac +1):ncol(A)]
  #facCov <- S[(nrow(S)-nfac + 1):nrow(S),(ncol(S)-nfac + 1):ncol(S)]
  #onemat <- matrix(1,nrow(facLoad),ncol(facLoad)) # matrix of ones for fac load deriv

  # from Cudeck 1993
  #h <- 0_00001
  # all possible parameters
  #pPos <- nvar*nfac + nfac*(nfac + 1)/2 + nvar*(nvar + 1)/2
  #pPos <- rep(1,length(par))

   grad_out <- rep(0,length(par))
  # h <- 0.00001
   h = sqrt(.Machine$double.eps)

   A_iter <- max(A)


  if(type=="none"){

    for(i in 1:length(par)){
      add <- rep(0,length(par))
      add[i] <- h
      ImpCovL = rcpp_RAMmult((par+add),A,S,S_fixed,A_fixed,A_est,S_est,F,I)[[1]]
      #ImpCovL = RAMmult((par+add),A,S,F,A_fixed,A_est,S_fixed,S_est)[[1]]
      ImpCovDot <- (ImpCovL - ImpCov)/h
      grad_out[i] <- 0.5 * trace(solve(ImpCov) %*% (ImpCov - SampCov) %*% solve(ImpCov) %*% ImpCovDot)
    }


  }  else if(type=="ridge"){ # not specified

    for(i in 1:length(par)){
      add <- rep(0,length(par))
      add[i] <- h
      ImpCovL = rcpp_RAMmult((par+add),A,S,S_fixed,A_fixed,A_est,S_est,F,I)[[1]]
      #ImpCovL = RAMmult((par+add),A,S,F,A_fixed,A_est,S_fixed,S_est)[[1]]
      ImpCovDot <- (ImpCovL - ImpCov)/h
      grad_out[i] <- 0.5 * (trace(solve(ImpCov) %*% (ImpCov - SampCov) %*% solve(ImpCov) %*% ImpCovDot)) +
        if(any(i==pars_pen)) 2*lambda*(max(Areg[A == i], Sreg[S==i])) else(0)
    }


  }  else if(type=="lasso"){
    #lasso

    for(i in 1:length(par)){
      add <- rep(0,length(par))
      add[i] <- h
      ImpCovL = rcpp_RAMmult((par+add),A,S,S_fixed,A_fixed,A_est,S_est,F,I)[[1]]
      #ImpCovL = RAMmult((par+add),A,S,F,A_fixed,A_est,S_fixed,S_est)[[1]]
      ImpCovDot <- (ImpCovL - ImpCov)/h
      grad_out[i] <- 0.5 * (trace(solve(ImpCov) %*% (ImpCov - SampCov) %*% solve(ImpCov) %*% ImpCovDot)) #+
      #  if(any(i==pars_pen)) lambda*sign(max(Areg[A == i], Sreg[S==i])) else(0)
    }


  }  else if(type=="enet"){ # not specified
    #elastic net

    for(i in 1:length(par)){
      add <- rep(0,length(par))
      add[i] <- 1 + h
      ImpCovL = RAMmult((par + add),A,S,F,A_fixed,A_est,S_fixed,S_est)[[1]]
      ImpCovDot <- (ImpCovL - ImpCov)/h
      grad_out[i] <- 0.5 * trace(solve(ImpCov) %*% (ImpCov - SampCov) %*% solve(ImpCov) %*% ImpCovDot)
    }


  }else if(type=="ols_lasso"){
    #lasso

    for(i in 1:length(par)){
      add <- rep(0,length(par))
      add[i] <- h
      ImpCovL = RAMmult((par + add),A,S,F,A_fixed,A_est,S_fixed,S_est)[[1]]
      ImpCovDot <- (ImpCovL - ImpCov)/h
      grad_out[i] <- 0.5 * trace(solve(ImpCov) %*% (ImpCov - SampCov) %*%
                                   solve(ImpCov) %*% ImpCovDot) + if(i <= A_iter) lambda*sign(Areg[A==i]) else(0)
    }


  }

  as.numeric(grad_out)


}