Consegui rodar o primeiro, instalando a mais nova versão do R no Windows. O segundo não consegui porque não forneceu o segundo arquivo 

dados2 <- read.table("HETE.csv", header = TRUE, sep = ";", dec = ",")

 

 

Prezados membros da lista,

Alguns dias postei um código não reproduzível sobre modelos para heterocedásticidade.

Peço desculpas e abaixo vai o código que imagino ser possível reproduzir.

Os dois principais problemas eram a instalação do pacote mcglm e o conjunto de dados.

O pacote pode ser instalado facilmente pelo github repository neste endereço

https://github.com/wbonat/mcglm
install_github("wbonat/mcglm", ref = "devel") # Saliento para instalar a versão devel !!

O conjunto de dados está anexado neste e-mail.

Além disso, o Luiz Leal estava com problemas pra ter acesso e pediu o .tar.gz que também vai em anexo.

Aproveito para salientar que o mcglm está em desenvolvimento, assim quaisquer dúvidas, criticas e/ou
sugestões serão muito bem vindas.

 

All the best!

# Heteroscedastic regression model -------------------------------------
# Author: Wagner Hugo Bonat LEG/UFPR -----------------------------------
# Date: 13/11/2016 -----------------------------------------------------

# Loading extra packages
install.packages("devtools")
require(devtools)

# Install mcglm package from github repository -------------------------
install_github("wbonat/mcglm", ref = "devel")
require(mcglm)


# Loading data set
Fenois <- c(337.311, 344.874, 342.353, 325.546, 333.950, 330.588, 328.067,
            328.067, 318.824, 331.429, 333.950, 334.790, 336.471, 338.151,
            342.353, 259.160, 252.437, 268.403, 265.882, 266.723, 287.731,
            88.571, 88.571,  90.252,  41.513,  52.437,  49.076,  88.571,
            88.571,  90.252,  64.202,  60.000,  61.681)
Cor <- factor(c(rep("ambar",6), rep("ambar_claro",3), rep("ambar",6),
                rep("ambar_claro",6),rep("branco",6),
                rep("extra_ambar_claro",3),rep("branco",3)))
dados <- data.frame("Fenois" = Fenois, "Cor" = Cor)
boxplot(dados$Fenois ~ dados$Cor)
tapply(dados$Fenois, dados$Cor, sd)

# Linear regression model- ---------------------------------------------
fit1 <- lm(Fenois ~ Cor, data = dados)
anova(fit1)
plot(residuals(fit1) ~ fitted(fit1))
plot(residuals(fit1) ~ Cor)

# Double Linear regression model ---------------------------------------
dados$id <- 1
fit2 <- mcglm(c(Fenois ~ Cor), list(mc_dglm(~ Cor, id = "id", data = dados)),
              covariance = "expm", data = dados)
summary(fit2)

summary(fit1)
summary(fit2)

cbind(coef(fit1), coef(fit2, type = "beta")$Estimates)
cbind(sqrt(diag(vcov(fit1))), coef(fit2, type = "beta", std.error = TRUE)$Std.error)

# Example 2 ------------------------------------------------------------
dados2 <- read.table("HETE.csv", header = TRUE, sep = ";", dec = ",")
with(dados2, boxplot(y ~ x))

# Note that, the A33 has no variance, so we need to remove this level.
dados2 <- dados2[which(dados2$x != "A33"),]
dados2$x <- droplevels(dados2$x)
tapply(dados2$y, dados2$x, sd)

# Linear regression model ----------------------------------------------
fit_lm <- lm(y ~ x, data = dados2)
summary(fit_lm)
plot(residuals(fit_lm) ~ fitted(fit_lm))
plot(residuals(fit_lm) ~ dados2$x)

# Double linear regression model ---------------------------------------
dados2$id <- 1
fit_dlm <- mcglm(c(y ~ x), list(mc_dglm(~ x, id = "id", data = dados2)),
              covariance = "expm", data = dados2)
summary(fit_dlm)

# Comparing estimates and standard errors ------------------------------
cbind(coef(fit_lm), coef(fit_dlm, type = "beta")$Estimates)
cbind(sqrt(diag(vcov(fit_lm))), coef(fit_dlm, type = "beta", std.error = TRUE)$Std.error)

 

fpr.br/r-br-guia) e forne�a c�digo m�nimo reproduz�vel.

_______________________________________________
R-br mailing list
R-br@listas.c3sl.ufpr.br
https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça código mínimo reproduzível.