R Under development (unstable) (2024-02-21 r85967 ucrt) -- "Unsuffered Consequences" Copyright (C) 2024 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library(lfe) Loading required package: Matrix > set.seed(65513) > options(lfe.threads=1,digits=4,warn=1) > x <- rnorm(52) > x2 <- rnorm(length(x)) > x3 <- 0.2*x + 0.1*x2 > ## create individual and firm > id <- factor(sample(30,length(x),prob=c(2,2,rep(1,28)),replace=TRUE)) > firm <- factor(sample(32,length(x),prob=c(2,2,rep(1,30)),replace=TRUE)) > year <- factor(sample(3,length(x),replace=TRUE)) > ## effects > id.eff <- rnorm(nlevels(id)) > firm.eff <- rnorm(nlevels(firm)) > year.eff <- rnorm(nlevels(year)) > ## left hand side > y <- x + 0.25*x2 + x3 + id.eff[id] + firm.eff[firm] + year.eff[year] + rnorm(length(x)) > > ## estimate and print result > est <- felm(y ~ x+x2+x3 |id+firm+year, exactDOF=TRUE) Warning in chol.default(mat, pivot = TRUE, tol = tol) : the matrix is either rank-deficient or not positive definite > > ## extract the group fixed effects > head(getfe(est)) effect obs comp fe idx id.1 4.4230 2 1 id 1 id.2 0.0000 5 1 id 2 id.3 0.0000 1 4 id 3 id.5 2.2393 2 1 id 5 id.6 -0.2621 3 1 id 6 id.7 3.5197 2 1 id 7 > tail(alpha <- getfe(est,ef='ln')) Warning in is.estimable(ef, obj$fe) : non-estimable function, largest error 0.7 in coordinate 3 ("id.3") Warning in getfe.kaczmarz(obj, se, ef = ef, bN = bN, robust = robust, cluster = cluster, : Supplied function seems non-estimable effect obs comp fe idx firm.29 1.8506 2 1 firm 29 firm.31 -2.9216 1 1 firm 31 firm.32 -1.2955 2 1 firm 32 year.1 -0.7471 20 5 year 1 year.2 -1.3103 16 5 year 2 year.3 0.1228 16 5 year 3 > # get the names to use below, just to make it easier > # lower precision in output > > nm <- rownames(alpha) > getfe(est,ef='zm',se=TRUE) effect obs comp se id.1 4.42299 2 1 1.0989 id.2 0.00000 5 1 0.0000 id.5 2.23926 2 1 0.7118 id.6 -0.26206 3 1 1.3249 id.7 3.51966 2 1 0.6274 id.10 0.54834 1 1 1.3315 id.11 4.52348 1 1 1.0906 id.12 -0.98009 1 1 0.8047 id.15 -0.64986 2 1 0.6681 id.16 1.70342 1 1 1.0450 id.17 -0.08610 1 1 0.9015 id.19 0.40255 3 1 0.4941 id.21 1.99314 5 1 0.3715 id.22 0.38724 2 1 0.7613 id.24 -1.21240 4 1 0.6319 id.27 -1.02366 2 1 1.2771 id.28 -0.75433 2 1 0.8389 firm.2 0.06196 1 1 0.9543 firm.4 3.11581 1 1 1.0748 firm.7 -0.16927 3 1 0.3753 firm.8 1.56558 5 1 0.8141 firm.11 -2.05081 3 1 0.4539 firm.13 -1.67562 3 1 0.6425 firm.15 -3.12297 3 1 0.7956 firm.16 -1.28607 3 1 0.6579 firm.17 -1.16917 4 1 0.4345 firm.18 1.57157 3 1 0.5064 firm.19 1.57518 2 1 0.4554 firm.20 3.98771 1 1 0.7692 firm.26 -0.24951 2 1 0.7042 firm.29 1.92129 2 1 0.6557 firm.31 -2.85090 1 1 0.8858 firm.32 -1.22477 2 1 0.7566 icpt.1 -1.75327 39 1 0.5419 id.14 0.09156 1 2 0.6734 id.30 0.00000 6 2 0.0000 firm.1 -1.88985 1 2 0.5527 firm.3 -0.36318 1 2 0.5458 firm.12 0.93403 1 2 0.5937 firm.22 -0.61133 2 2 0.5910 firm.24 0.52948 1 2 0.5169 firm.25 1.40085 1 2 0.6002 icpt.2 -0.89466 7 2 0.4346 id.13 0.79325 1 3 1.0520 id.29 0.00000 4 3 0.0000 firm.9 -1.35501 2 3 0.3693 firm.14 -0.31142 1 3 0.4128 firm.27 1.66643 2 3 0.2922 icpt.3 -0.84079 5 3 0.3207 id.3 0.00000 1 4 0.0000 firm.23 0.00000 1 4 0.0000 icpt.4 -2.77636 1 4 0.6033 year.1 0.00000 20 5 0.0000 year.2 -0.56318 16 5 0.4197 year.3 0.86994 16 5 0.7280 > getfe(est,ef='zm2',se=TRUE) effect obs comp se id.1 3.55407 2 1 1.0564 id.2 -0.86892 5 1 0.3668 id.5 1.37034 2 1 0.5317 id.6 -1.13097 3 1 0.9327 id.7 2.65074 2 1 0.5574 id.10 -0.32058 1 1 0.9228 id.11 3.65457 1 1 0.7661 id.12 -1.84900 1 1 0.8004 id.15 -1.51878 2 1 0.6700 id.16 0.83450 1 1 1.1504 id.17 -0.95502 1 1 0.7420 id.19 -0.46637 3 1 0.4815 id.21 1.12422 5 1 0.3946 id.22 -0.48168 2 1 0.9089 id.24 -2.08132 4 1 0.4684 id.27 -1.89258 2 1 0.9914 id.28 -1.62324 2 1 0.4701 firm.2 0.06196 1 1 0.8749 firm.4 3.11581 1 1 0.9183 firm.7 -0.16927 3 1 0.3346 firm.8 1.56558 5 1 0.7551 firm.11 -2.05081 3 1 0.5050 firm.13 -1.67562 3 1 0.5538 firm.15 -3.12297 3 1 0.8583 firm.16 -1.28607 3 1 0.6305 firm.17 -1.16917 4 1 0.4343 firm.18 1.57157 3 1 0.5186 firm.19 1.57518 2 1 0.5217 firm.20 3.98771 1 1 0.7220 firm.26 -0.24951 2 1 0.6608 firm.29 1.92129 2 1 0.5953 firm.31 -2.85090 1 1 0.9308 firm.32 -1.22477 2 1 0.7421 icpt.1 -0.88436 39 1 0.3651 id.14 0.04578 1 2 0.3438 id.30 -0.04578 6 2 0.3438 firm.1 -1.88985 1 2 0.5069 firm.3 -0.36318 1 2 0.5806 firm.12 0.93403 1 2 0.5308 firm.22 -0.61133 2 2 0.5728 firm.24 0.52948 1 2 0.5473 firm.25 1.40085 1 2 0.6711 icpt.2 -0.84888 7 2 0.4918 id.13 0.39662 1 3 0.4530 id.29 -0.39662 4 3 0.4530 firm.9 -1.35501 2 3 0.3890 firm.14 -0.31142 1 3 0.4055 firm.27 1.66643 2 3 0.2644 icpt.3 -0.44417 5 3 0.4349 id.3 0.00000 1 4 0.0000 firm.23 0.00000 1 4 0.0000 icpt.4 -2.77636 1 4 0.6537 year.1 0.00000 20 5 0.0000 year.2 -0.56318 16 5 0.3755 year.3 0.86994 16 5 0.6874 > ef <- function(v,addnames) { + names(v) <- nm + w <- c(v['id.2']-v['id.1'],exp(v['id.2']-v['id.1']), + v['id.2']+v['firm.2']+v['year.1'],v['id.5']+v['firm.4']+v['year.2']) + if(addnames) names(w) <-c('id2-id1','exp(id2-id1)','id2+f2+y1','id5+f4+y2') + w + } > getfe(est,ef=ef,se=TRUE) effect se id2-id1 -4.423 1.037558 exp(id2-id1) 0.012 0.002014 id2+f2+y1 -1.691 1.266594 id5+f4+y2 3.039 1.556251 > > # test whether we have estimable functions > > R <- est$r.residuals - est$residuals > > cat('myef :',is.estimable(ef,est$fe,R),'\n') myef : TRUE > for(n in c('ref','zm','zm2','ln')) { + cat(n,':',is.estimable(efactory(est,n),est$fe,R),'\n') + } ref : TRUE zm : TRUE zm2 : TRUE Warning in is.estimable(efactory(est, n), est$fe, R) : non-estimable function, largest error 0.6 in coordinate 41 ("firm.23") ln : FALSE > > proc.time() user system elapsed 1.67 0.15 1.82