Это будет длинный вопрос с довольно большими частями кода, но, пожалуйста, помогите. Вы можете пропустить большую часть этого, я просто должен был положить все из-за характера моей проблемы, так что извините за это. Это длинное чтение, но также быстрое.
Я пытаюсь оценить модель SUR в R с помощью пакета systemfit. Это довольно простая задача. Я загружаю данные .xlsx в R, определяю систему уравнений как список объектов «формулы класса», устанавливаю method = "SUR", и это работает. Вот мой код, который делает это без ошибок.
#FUNCTIONS THAT I NEED TO TRANSFORM MY DATA BEFORE FITTING THE MODEL
testdata <- read_excel("testdata.xlsx")
attach(testdata)
differenceData <- function(var, l, order){
#function that returns difference of a variable for a given lag and order
zeros <- c(rep(0, l))
if (order == 0) {
return (var)
} else if (order == 1) {
return(c(zeros,diff(var,lag = l, differences = order)))
} else if (order == 2) {
return(c(zeros,diff(var,lag = l, differences = order)))
}
}
lagSeries <- function(var, order){
#function that takes a dataset variable and lags it
if (order == 1) {
zeros = 0
var <- var[-order]
var <- c(var, zeros)
return (var)
} else {
zeros = rep(0, order)
for (i in 1:order) {
var <- var[-1]
}
var <- c(var, zeros)
return(var)
}
}
#Extracting variables necessary for SUR model
Xset <- cbind(gdp_index_sa, consumer_price_index, active_interest_rate, employment_inedx)
Y <- differenceData(npl_consumer, 1, 1)
#Creating equations for my model
main <- Y ~ Xset #first equation in SUR model, multiple regression
ar1 <- Xset[,1] ~ lagSeries(Xset[,1], 1) #2nd equations - ar(1) quation
ar2 <- Xset[,2] ~ lagSeries(Xset[,2], 1) #3rd equation - ar(1) equation
ar3 <- Xset[,3] ~ lagSeries(Xset[,3], 1) #4th equation - ar(1) equation
ar4 <- Xset[,4] ~ lagSeries(Xset[,4], 1) #5th quation - ar(1) equation
system <- list(main,ar1,ar2,ar3,ar4) #list of equations for estimation, a 'system' argument for #systemfit
model <- systemfit :: systemfit(system, method = "SUR") #fitting the model
summary(model)
Это прекрасно работает. Вот результат (не стесняйтесь, это хорошо)
SUR estimates for 'eq1' (equation 1)
Model Formula: Y ~ Xset
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.034084659 0.014006080 2.43356 0.0158478 *
Xsetgdp_index_sa 0.000773736 0.001652541 0.46821 0.6401546
Xsetconsumer_price_index -0.011804724 0.004075654 -2.89640 0.0042036 **
Xsetactive_interest_rate -0.000500245 0.001624472 -0.30794 0.7584526
Xsetemployment_inedx -0.001988833 0.001566473 -1.26963 0.2057236
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.04699 on 196 degrees of freedom
Number of observations: 201 Degrees of Freedom: 196
SSR: 0.432782 MSE: 0.002208 Root MSE: 0.04699
Multiple R-Squared: 0.064292 Adjusted R-Squared: 0.045196
SUR estimates for 'eq2' (equation 2)
Model Formula: Xset[, 1] ~ lagSeries(Xset[, 1], 1)
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.574679 0.282222 5.57958 7.8112e-08 ***
lagSeries(Xset[, 1], 1) 0.708644 0.046616 15.20175 < 2.22e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.927956 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 739.685556 MSE: 3.717013 Root MSE: 1.927956
Multiple R-Squared: 0.545061 Adjusted R-Squared: 0.542775
SUR estimates for 'eq3' (equation 3)
Model Formula: Xset[, 2] ~ lagSeries(Xset[, 2], 1)
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.1482021 0.0508081 2.9169 0.0039415 **
lagSeries(Xset[, 2], 1) 0.9345933 0.0217035 43.0618 < 2.22e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.287545 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 16.453722 MSE: 0.082682 Root MSE: 0.287545
Multiple R-Squared: 0.902375 Adjusted R-Squared: 0.901885
SUR estimates for 'eq4' (equation 4)
Model Formula: Xset[, 3] ~ lagSeries(Xset[, 3], 1)
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.0615835 0.0606266 -1.01578 0.31097
lagSeries(Xset[, 3], 1) 0.9297783 0.0279298 33.28983 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.829464 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 136.914024 MSE: 0.68801 Root MSE: 0.829464
Multiple R-Squared: 0.843468 Adjusted R-Squared: 0.842682
SUR estimates for 'eq5' (equation 5)
Model Formula: Xset[, 4] ~ lagSeries(Xset[, 4], 1)
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.0144565 0.0558382 0.2589 0.79598
lagSeries(Xset[, 4], 1) 0.9918317 0.0129945 76.3273 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.529632 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 55.821572 MSE: 0.28051 Root MSE: 0.529632
Multiple R-Squared: 0.96666 Adjusted R-Squared: 0.966492
Вы видите коэффициенты, R Square et c - все нормальные значения. Проблема возникает, когда я хочу создать этот «системный» аргумент другим способом. Используя функцию. Моя идея - написать функцию, которая принимает два аргумента.
Y - зависимая переменная, как я ее здесь называю
Xset - набор независимых переменных
, так что вместо
system <- list(main,ar1,ar2,ar3,ar4)
я могу сделать что-то вроде
system <- define_system_of_equations(Y, Xset)
Вот тут у меня проблемы. Вот мой код, который делает это
getSystem2 <- function(Y, Xset){
#function that returns a list that should be a 'system' argument for systemfit
main <- Y ~ Xset
system <- list(main)
for (i in 1:ncol(Xset)) {
Xt <- Xset[,i] #I tries withoud as.numeric(unlist()) and I got same results
Xt1 <- lagSeries(Xt, 1)
ar1 <- Xt ~ Xt1
system[[i+1]] <- ar1
}
return(system)
}
system_app <- getSystem2(dep, Xset)
model_app <- systemfit::systemfit(system_app, method = "SUR")
summary(model_app)
Это не показывает ошибок, но результаты в сводке модели просто неверны. Я получаю отрицательные значения R в квадрате для моих уравнений ar (1), и некоторые из моих данных отсутствуют (выводит NA). Вот результат (проверьте нижнюю часть и сравните ее с нижней частью предыдущего вывода)
SUR estimates for 'eq1' (equation 1)
Model Formula: Y ~ Xset
<environment: 0x210d4520>
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.021438350 0.014036301 1.52735 0.128286
Xsetgdp_index_sa 0.001483992 0.001658170 0.89496 0.371908
Xsetconsumer_price_index -0.010368511 0.004082449 -2.53978 0.011868 *
Xsetactive_interest_rate -0.000505426 0.001630243 -0.31003 0.756867
Xsetemployment_inedx -0.002299415 0.001569858 -1.46473 0.144597
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.047456 on 196 degrees of freedom
Number of observations: 201 Degrees of Freedom: 196
SSR: 0.441411 MSE: 0.002252 Root MSE: 0.047456
Multiple R-Squared: 0.045637 Adjusted R-Squared: 0.02616
SUR estimates for 'eq2' (equation 2)
Model Formula: Xt ~ Xt1
<environment: 0x210d4520>
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.83430e-01 4.28569e-19 2.29468e+18 < 2.22e-16 ***
Xt1 9.83467e-01 5.17256e-35 1.90132e+34 < 2.22e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.085128 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 234.323071 MSE: 1.177503 Root MSE: 1.085128
Multiple R-Squared: 0.860048 Adjusted R-Squared: 0.859344
SUR estimates for 'eq3' (equation 3)
Model Formula: Xt ~ Xt1
<environment: 0x210d4520>
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.83430e-01 5.14996e-19 1.90959e+18 < 2.22e-16 ***
Xt1 9.83467e-01 NA NA NA
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.085128 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 234.323071 MSE: 1.177503 Root MSE: 1.085128
Multiple R-Squared: 0.860048 Adjusted R-Squared: 0.859344
SUR estimates for 'eq4' (equation 4)
Model Formula: Xt ~ Xt1
<environment: 0x210d4520>
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.25196e-14 NA NA NA
Xt1 3.55476e-15 3.33565e-19 10656.9 < 2.22e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 4.317739 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 3709.931822 MSE: 18.642873 Root MSE: 4.317739
Multiple R-Squared: -1.215801 Adjusted R-Squared: -1.226936
SUR estimates for 'eq5' (equation 5)
Model Formula: Xt ~ Xt1
<environment: 0x210d4520>
Estimate Std. Error t value Pr(>|t|)
(Intercept) -5.49310e-16 NA NA NA
Xt1 1.55969e-16 NA NA NA
Residual standard error: 4.317739 on 199 degrees of freedom
Number of observations: 201 Degrees of Freedom: 199
SSR: 3709.931822 MSE: 18.642873 Root MSE: 4.317739
Multiple R-Squared: -1.215801 Adjusted R-Squared: -1.226936
Обратите внимание, как R-квадрат внезапно становится отрицательным и как у меня есть некоторые пропущенные значения для коэф. значения et c. Это действительно странно, учитывая, что одни и те же данные используются в обеих аппроксимациях, и ничего не показывает ошибок.
Причина, по которой мне нужно это сделать с помощью функции, заключается в том, что я создаю блестящее приложение, в котором пользователь выбирает переменные для модели SUR и затем результаты форматируются и распечатываются для него в приложении. С этой целью все правильно, форматирование работает в приложении (у меня есть другие функции, чтобы делать все это), но эти цифры просто неверны.
Примечание: это не код приложения. Использование getSystem2 () просто имитирует то, что сделано в приложении, и я получаю те же проблемы локально.
Почему это происходит? Любые идеи? Я застрял на этом некоторое время, и мне действительно нужна помощь. Просто не могу понять это. Все отлично работает без функции, но это просто не для меня.
Я полагаю, что это воспроизводимо с любыми имеющимися у вас данными.
Есть идеи? Спасибо.