R: избегая summary.plm

Я использую R для запуска симуляции Монте-Карло, изучая производительность оценщиков панельных данных. Поскольку я буду запускать большое количество испытаний, мне нужно получить по крайней мере достойную производительность от моего кода.

Использование Rprof в 10 испытаниях моей симуляции показывает, что значительная часть времени тратится на вызовы summary.plm. Первые несколько строк Rprofsummary представлены ниже:

$by.total
                            total.time total.pct self.time self.pct
"trial"                          54.48     100.0      0.00      0.0
"coefs"                          53.90      98.9      0.06      0.1
"model.matrix"                   36.72      67.4      0.10      0.2
"model.matrix.pFormula"          35.98      66.0      0.06      0.1
"summary"                        33.82      62.1      0.00      0.0
"summary.plm"                    33.80      62.0      0.08      0.1
"r.squared"                      29.00      53.2      0.02      0.0
"FUN"                            24.84      45.6      7.52     13.8

Я вызываю summary в своем коде, потому что мне нужно получить стандартные ошибки оценок коэффициентов, а также сами коэффициенты (которые я мог бы получить только из объекта plm). Мой звонок выглядит как

regression <- plm(g ~ y0 + Xit, data=panel_data, model=model, index=c("country","period"))

coefficients_estimated <- summary(regression)$coefficients[,"Estimate"]
ses_estimated <- summary(regression)$coefficients[,"Std. Error"]

У меня есть непреодолимое чувство, что это огромная трата процессорного времени, но я недостаточно знаю, как R делает что-то, чтобы избежать вызова сводки. Я был бы признателен за любую информацию о том, что здесь происходит за кулисами, или о каком-либо способе сокращения времени, необходимого для выполнения этого.


person Wilduck    schedule 11.04.2011    source источник
comment
Хммм, первая мысль - объединить два сводных вызова. Я собираюсь попробовать это, и посмотреть, поможет ли это...   -  person Wilduck    schedule 11.04.2011
comment
Объединение двух сводных вызовов действительно немного помогает, но summary.plm по-прежнему остается серьезной проблемой.   -  person Wilduck    schedule 11.04.2011


Ответы (3)


Вам просто нужно заглянуть внутрь plm:::summary.plm, чтобы увидеть, что он делает. Когда вы это сделаете, вы увидите, что ваши две строки, вызывающие summary() для вашей модели, могут быть заменены на:

coefficients_estimated <- coef(regression)
ses_estimated <- sqrt(diag(vcov(regression)))

Например:

require(plm)
data("Produc", package = "plm")
zz <- plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp, 
          data = Produc, index = c("state","year"))

summary(zz) дает:

> summary(zz)
Oneway (individual) effect Within Model

....

Coefficients :
             Estimate  Std. Error t-value  Pr(>|t|)    
log(pcap) -0.02614965  0.02900158 -0.9017    0.3675    
log(pc)    0.29200693  0.02511967 11.6246 < 2.2e-16 ***
log(emp)   0.76815947  0.03009174 25.5273 < 2.2e-16 ***
unemp     -0.00529774  0.00098873 -5.3582 1.114e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
....

и две строки, которые я показал, возвращаются для zz:

> coef(zz)
   log(pcap)      log(pc)     log(emp)        unemp 
-0.026149654  0.292006925  0.768159473 -0.005297741 
> sqrt(diag(vcov(zz)))
   log(pcap)      log(pc)     log(emp)        unemp 
0.0290015755 0.0251196728 0.0300917394 0.0009887257

На самом деле вы не предоставляете достаточно информации (например, ваш код симуляции или полный вывод из Rprof()), чтобы сказать, поможет ли это - определенно не похоже, что в summary() тратится огромное количество времени; FUN намного дороже, чем все остальное, что вы показываете, а из элементов, которые вы показываете, r.squared() — единственный, который появляется в plm:::summary.plm(), и, похоже, это не занимает вообще никакого времени.

Итак, еще неизвестно, ускорит ли вышеизложенное процесс.

person Gavin Simpson    schedule 11.04.2011
comment
Вы правы, это незначительно влияет на скорость. Как я отметил в своем комментарии к вопросу, у меня было значительное увеличение от сокращения до одного звонка до сводки (сокращение времени примерно на 30%). Я не подумал посмотреть, что делает summary.plm, поэтому я ценю это понимание для будущих проблем. Что касается вашего совета посмотреть на FUN, я не уверен, что FUN имеет в виду. Я бы подумал, что аргумент применить, но у меня нет циклов применения ... Есть ли способ легко узнать, что такое FUN? В любом случае ответ принят. Спасибо. - person Wilduck; 13.04.2011
comment
@Wilduck FUN обычно является именем аргумента функции, которая должна применяться к фрагментам данных. У нас нет кода вашей симуляции, поэтому трудно сказать, что такое FUN, но используете ли вы какие-либо из семейств apply(), sapply(), lapply() и т. д.? - person Gavin Simpson; 13.04.2011
comment
да, примерно так я и думал. Я не использую какое-либо семейство приложений для больших блоков данных, поэтому это должно быть что-то другое. Я не собираюсь заставлять вас просматривать код моей симуляции до тех пор, пока я не пройдусь по нему гребешком с мелкими зубьями. Спасибо за вашу помощь. - person Wilduck; 14.04.2011

Если вы хотите пойти дальше, взгляните на фактический код функции plm:::plm. Вы заметите, что перед окончательным вызовом plm:::plm.fit выполняется много проверок аргументов. Вы можете (если действительно хотите) сразу перейти к plm.fit.

Один последний момент. Вы упомянули, что ваша проблема — это симуляция Монте-Карло. Можете ли вы использовать параллельные вычисления для увеличения скорости?

person csgillespie    schedule 11.04.2011
comment
Хорошо, жаль, что эта функция не связана с ?plm, поскольку мне было интересно, есть ли мощное вычислительное ядро, которое @Wilduck мог бы использовать. Я бросил беглый взгляд на справку по plm(). - person Gavin Simpson; 11.04.2011
comment
Как вы заметили, я думаю, что просмотр FUN может быть лучшим выбором для @wilduck. - person csgillespie; 11.04.2011
comment
Я планирую использовать параллельные вычисления, поскольку монте-карло смущающе параллелен, я думал, что foreach будет работать довольно хорошо. - person Wilduck; 13.04.2011

Просто используйте coeftest(zz). coeftest находится в пакете lmtest; это даст вам коэффициенты и стандартные ошибки из plm объектов намного быстрее, чем summary.plm.

person user697473    schedule 04.03.2012