Re: [R-br] Estatística circular - gráfico
Segue: Dados AMBIENTE mess angle DISTANCIA ABUND RIQUEZA NASC 1 0 0 2 2 NASC 2 30 0 3 1 NASC 3 60 0 4 3 NASC 4 90 0 8 5 NASC 5 120 0 0 0 NASC 6 150 0 0 0 NASC 7 180 0 31 1 NASC 8 210 0 10 3 NASC 9 240 0 18 3 NASC 10 270 0 7 5 NASC 11 300 0 2 2 NASC 12 330 0 1 1 ############################################################### ################### início da análise ############################## ############################################################### library(circular) nascente=read.table("nasc.csv",h=T) nascente nascente0=nascente[which(nascente$DISTANCIA==0),] nascente0 meses<-c("Mai","Jun","Jul","Ago","Set","Out","Nov","Dez","Jan","Fev","Mar","Abr") #convertendo os dados para circular, usando radianos nasc.c0<-circular(rep(rad(nascente0$angle),nascente0$ABUND),units=c("radians"),rotation="clock") #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Angulo médio + Desvio padrão angular%%%%%%%%%%%%%%%%%%%% # %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% media=mean(nasc.c0) sqrt(-2*log(rho.circular(nasc.c0))) #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% #### Gráfico ############ par(mfrow=c(1,3)) par(mar=c(0,0,2,0)) plot(nasc.c0,axes=F,main="0") axis.circular(at=circular(rad(nascente0$angle)),labels=meses,units=c("radians"), rotation="clock") rose.diag(nasc.c0,add=T,axes=F) arrows.circular(median(nasc.c0),y=rho.circular(nasc.c0),lwd=4,lty=1,angle=90,len=0.0,col="green") #aqui eu inseri a mediana arrows.circular(mean(nasc.c0),y=rho.circular(nasc.c0),lwd=4,lty=1,angle=90,len=0.0,col="red") #aqui eu inseri a média #agora preciso inserir os desvios, ou erros padrão. ########################### FIM ###################################### Obrigado Em 8 de junho de 2016 12:00, <r-br-request@listas.c3sl.ufpr.br> escreveu:
Enviar submissões para a lista de discussão R-br para r-br@listas.c3sl.ufpr.br
Para se cadastrar ou descadastrar via WWW, visite o endereço https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br ou, via email, envie uma mensagem com a palavra 'help' no assunto ou corpo da mensagem para r-br-request@listas.c3sl.ufpr.br
Você poderá entrar em contato com a pessoa que gerencia a lista pelo endereço r-br-owner@listas.c3sl.ufpr.br
Quando responder, por favor edite sua linha Assunto assim ela será mais específica que "Re: Contents of R-br digest..."
Tópicos de Hoje:
1. Re: Estatística circular - gráfico (Éder Comunello) 2. Erro modelo - glm (Andre Oliveira) 3. Re: Erro modelo - glm (Marcus Nunes) 4. Re: Erro modelo - glm (Leonardo Ferreira Fontenelle)
----------------------------------------------------------------------
Message: 1 Date: Tue, 7 Jun 2016 14:39:28 -0400 From: Éder Comunello <comunello.eder@gmail.com> To: "r-br@listas.c3sl.ufpr.br" <r-br@listas.c3sl.ufpr.br> Subject: Re: [R-br] Estatística circular - gráfico Message-ID: <CABmC8gn-f8sdoHs0-Hc=zAGTmd=iYjMHizFe= X_HK6JYa4Vobw@mail.gmail.com> Content-Type: text/plain; charset="utf-8"
Cassiano, boa tarde!
Talvez seja o caso de postar seu código com um dataset de exemplo.
================================================ Éder Comunello Researcher at Brazilian Agricultural Research Corporation (Embrapa) DSc in Agricultural Systems Engineering (USP/Esalq) MSc in Environ. Sciences (UEM), Agronomist (UEM) --- Embrapa Agropecuária Oeste, Dourados, MS, Brazil |<O>| ================================================ GEO, -22.2752, -54.8182, 408m UTC-04:00 / DST: UTC-03:00
Em 7 de junho de 2016 10:11, Cassiano <cassianosr@gmail.com> escreveu:
Obrigado pela ajuda. Realmente o ggplot2 produz gráficos mais interessantes. Mas ainda não achei a opção para plotar a média e os desvios angulares. Se puderem me dar mais instruções, ficarei agradecido.
Abraços Cassiano
-------------- Próxima Parte ---------- Um anexo em HTML foi limpo... URL: < http://listas.inf.ufpr.br/pipermail/r-br/attachments/20160607/2bc39a6a/attac...
------------------------------
Message: 2 Date: Tue, 7 Jun 2016 23:26:58 +0000 (UTC) From: Andre Oliveira <andreolsouza@yahoo.com.br> To: R-br Lista <r-br@listas.c3sl.ufpr.br> Subject: [R-br] Erro modelo - glm Message-ID: <1694617424.92644.1465342018577.JavaMail.yahoo@mail.yahoo.com> Content-Type: text/plain; charset="utf-8"
boa noite, alguém poderia ajudar descobrir onde está o erro?
remove(list=ls())diurno=c(68,36,40,30) noturno=c(46,26,55,22) quadra= c("A", "B", "G","D") dados=data.frame(diurno,noturno,quadra) fit = glm(cbind(diurno, noturno) ~ factor(quadra), family = poisson, data=dados) Error in x[good, , drop = FALSE] : (subscript) subscrito lógico muito longo obrigado André Oliveira Souza. Graduação em Matemática, mestrado em estatística aplicada.Instituto Federal de Educação, Ciência e Tecnologia do Espirito Santo. IFES -------------- Próxima Parte ---------- Um anexo em HTML foi limpo... URL: < http://listas.inf.ufpr.br/pipermail/r-br/attachments/20160607/7d2e5cef/attac...
------------------------------
Message: 3 Date: Tue, 7 Jun 2016 22:37:17 -0300 From: Marcus Nunes <marcus.nunes@gmail.com> To: r-br@listas.c3sl.ufpr.br, Andre Oliveira <andreolsouza@yahoo.com.br> Subject: Re: [R-br] Erro modelo - glm Message-ID: < CA+QGQvvGr8Rg5w+mO4tNmRiG1du+nF0Z83XCZUrDfYyhgZdTww@mail.gmail.com> Content-Type: text/plain; charset="utf-8"
Ajuste os modelos individualmente:
diurno=c(68,36,40,30) noturno=c(46,26,55,22) quadra= c("A", "B", "G","D") dados=data.frame(diurno,noturno,quadra) fit1 = glm(diurno ~ factor(quadra), family = poisson, data=dados) fit2 = glm(noturno ~ factor(quadra), family = poisson, data=dados)
Isto deve resolver o teu problema.
2016-06-07 20:26 GMT-03:00 Andre Oliveira <andreolsouza@yahoo.com.br>:
boa noite, alguém poderia ajudar descobrir onde está o erro?
remove(list=ls()) diurno=c(68,36,40,30) noturno=c(46,26,55,22) quadra= c("A", "B", "G","D") dados=data.frame(diurno,noturno,quadra) fit = glm(cbind(diurno, noturno) ~ factor(quadra), family = poisson, data=dados)
*Error in x[good, , drop = FALSE] : * * (subscript) subscrito lógico muito longo*
obrigado
André Oliveira Souza. Graduação em Matemática, mestrado em estatística aplicada.Instituto Federal de Educação, Ciência e Tecnologia do Espirito Santo. IFES
_______________________________________________ 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.
-- Marcus Nunes http://marcusnunes.me/ -------------- Próxima Parte ---------- Um anexo em HTML foi limpo... URL: < http://listas.inf.ufpr.br/pipermail/r-br/attachments/20160607/a1cd75b1/attac...
------------------------------
Message: 4 Date: Wed, 08 Jun 2016 11:16:04 -0300 From: Leonardo Ferreira Fontenelle <leonardof@leonardof.med.br> To: r-br@listas.c3sl.ufpr.br Subject: Re: [R-br] Erro modelo - glm Message-ID: <1465395364.4116631.631632273.0B1F3BA2@webmail.messagingengine.com
Content-Type: text/plain; charset="utf-8"
boa noite, alguém poderia ajudar descobrir onde está o erro?
remove(list=ls()) diurno=c(68,36,40,30) noturno=c(46,26,55,22) quadra= c("A", "B", "G","D") dados=data.frame(diurno,noturno,quadra) fit = glm(cbind(diurno, noturno) ~ factor(quadra), family = poisson, data=dados) *Error in x[good, , drop = FALSE] :
Em Ter 7 jun. 2016, às 20:26, Andre Oliveira escreveu: *
* (subscript) subscrito lógico muito longo*
Andre, conforme está escrito na documentação do glm, a resposta pode ser fornecida como duas colunas (sucessos e fracassos) *apenas para regressões binomiais e quasi-binomiais*; não para Poisson.
Supondo que seus dados sejam o resumo de 68 + 46 = 114 observações para a quadra A, 36 + 26 = 62 observações para a quadra B, e por aí em diante; e que você queira calcular uma razão de prevalência / risco / outra proporção entre as quadras, o jeito mais óbvio de fornecer seus dados para a função é em formato longo:
dados <- data.frame(horariodiurno = c(1, 0, 1, 0, 1, 0, 1, 0), quadra = c("A", "A", "B", "B", "C", "C", "D", "D")) dados <- dados[rep(1:nrow(dados), c(68, 46, 36, 26, 40, 55, 30, 22)), ] fit = glm(horariodiurno ~ quadra, family = quasipoisson, data = dados)
Repare também que, se você quer trabalhar com razão de proporções, a distribuição Poisson em princípio não se aplica, daí a necessidade de você especificar uma quasipoisson ou então um estimador de variância consistente com heterocedasticidade (pacote sandwich).
Por favor, use o código com cuidado porque não estou com meu computador aqui para testá-lo.
Atenciosamente,
Leonardo Ferreira Fontenelle[1]
Links:
1. http://lattes.cnpq.br/9234772336296638 -------------- Próxima Parte ---------- Um anexo em HTML foi limpo... URL: < http://listas.inf.ufpr.br/pipermail/r-br/attachments/20160608/7af199f3/attac...
------------------------------
Subject: Legenda do Digest
_______________________________________________ R-br mailing list R-br@listas.c3sl.ufpr.br https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
------------------------------
Fim da Digest R-br, volume 66, assunto 7 ****************************************
-- ======================================= Cassiano S. Rosa CRBio: 104117/04-D Universidade Federal do Triângulo Mineiro (UFTM), *campus* Iturama-MG Avenida Rio Paranaíba, 1241, Centro, 38280-000 Tel: (34)3415-2512 Lattes <http://lattes.cnpq.br/3030825329622014> Scholar Google <http://scholar.google.com/citations?user=JmELzOMAAAAJ> Academia.edu <https://uftm.academia.edu/CassianoRosa> Researchgate <https://www.researchgate.net/profile/Cassiano_Rosa> =======================================
Cassiano, bom dia! Imagino que o que você esteja tentando adicionar ao gráfico é uma linha da média + erro tal qual essa que aparece no gráfico abaixo: Infelizmente, não encontrei como reproduzir esse efeito no R. O melhor que consegui foi marcar os limites com novas linhas. [image: Imagem inline 1] ### <code r> library(circular) nascente <- read.table(text= "AMBIENTE mes angle DISTANCIA ABUND RIQUEZA NASC 1 0 0 2 2 NASC 2 30 0 3 1 NASC 3 60 0 4 3 NASC 4 90 0 8 5 NASC 5 120 0 0 0 NASC 6 150 0 0 0 NASC 7 180 0 31 1 NASC 8 210 0 10 3 NASC 9 240 0 18 3 NASC 10 270 0 7 5 NASC 11 300 0 2 2 NASC 12 330 0 1 1", head=T) nascente nascente0=nascente[which(nascente$DISTANCIA==0),]; nascente0 meses <- c("Mai","Jun","Jul","Ago","Set","Out","Nov","Dez","Jan","Fev","Mar","Abr") #convertendo os dados para circular, usando radianos nasc.c0 <- circular(rep(rad(nascente0$angle),nascente0$ABUND), units=c("radians"), rota="clock") nasc.c0 plot(nasc.c0, stack=T) # Angulo médio + Desvio padrão angular media <- mean(nasc.c0) dp <- sd(nasc.c0) n <- length(nasc.c0) err <- dp/sqrt(n) len <- rho.circular(nasc.c0) cbind(media, dp, n, err, len) sd(nasc.c0); sqrt(-2*log(rho.circular(nasc.c0))) # Gráfico # par(mfrow=c(1,3)) # par(mar=c(0,0,2,0)) plot(nasc.c0,axes=F,main=NA) axis.circular(at=circular(rad(nascente0$angle)), labels=meses, units=c("radians"), rota="clock") rose.diag(nasc.c0, add=T, axes=F, col=8, bor=8, lwd=3) # arrows.circular(median(nasc.c0),y=rho.circular(nasc.c0),lwd=3,lty=1,angle=30,len=0.1,col="green") #aqui eu inseri a mediana arrows.circular(mean(nasc.c0),y=len, lwd=2, lty=1, angle=90, len=0.1, col=2) arrows.circular(media+err, y=1, lwd=1, lty=3, angle=0, len=0.1, col=2) arrows.circular(media-err, y=1, lwd=1, lty=3, angle=0, len=0.1, col=2) ### </code> ================================================ Éder Comunello Researcher at Brazilian Agricultural Research Corporation (Embrapa) DSc in Agricultural Systems Engineering (USP/Esalq) MSc in Environ. Sciences (UEM), Agronomist (UEM) --- Embrapa Agropecuária Oeste, Dourados, MS, Brazil |<O>| ================================================ GEO, -22.2752, -54.8182, 408m UTC-04:00 / DST: UTC-03:00
Cassiano, olá novamente! Atualizei o código, pois percebi que havia um erro na atribuição dos labels (meses) que estavam "counter-clock". [image: Imagem inline 1] ### <code r> library(circular) nascente <- read.table(text= "AMBIENTE mes angle DISTANCIA ABUND RIQUEZA NASC 1 0 0 2 2 NASC 2 30 0 3 1 NASC 3 60 0 4 3 NASC 4 90 0 8 5 NASC 5 120 0 0 0 NASC 6 150 0 0 0 NASC 7 180 0 31 1 NASC 8 210 0 10 3 NASC 9 240 0 18 3 NASC 10 270 0 7 5 NASC 11 300 0 2 2 NASC 12 330 0 1 1", head=T) nascente$meses <- c("Mai","Jun","Jul","Ago","Set","Out","Nov","Dez","Jan","Fev","Mar","Abr") nascente #convertendo os dados para circular, usando radianos nasc.c0 <- circular(rep(rad(nascente$angle),nascente$ABUND), units=c("radians"), rota="clock") nasc.c0 plot(nasc.c0, stack=T) deg(nasc.c0) # Angulo médio + Desvio padrão angular media <- mean(nasc.c0) dp <- sd(nasc.c0) n <- length(nasc.c0) err <- dp/sqrt(n) len <- rho.circular(nasc.c0) cbind(media, dp, n, err, len) sd(nasc.c0); sqrt(-2*log(rho.circular(nasc.c0))) # Gráfico # par(mfrow=c(1,3)) # par(mar=c(0,0,2,0)) plot(nasc.c0,axes=F, main=NA, stack=F, col=3) axis.circular(at=circular(rad(nascente$angle), rota="clock"), labels=nascente$meses, units=c("radians"), col=4) rose.diag(nasc.c0, add=T, axes=F, col=8, bor=8, lwd=3) # arrows.circular(median(nasc.c0),y=rho.circular(nasc.c0),lwd=3,lty=1,angle=30,len=0.1,col="green") #aqui eu inseri a mediana arrows.circular(mean(nasc.c0),y=len, lwd=2, lty=1, angle=90, len=0.1, col=2) arrows.circular(media+err, y=1, lwd=1, lty=3, angle=0, len=0.1, col=2) arrows.circular(media-err, y=1, lwd=1, lty=3, angle=0, len=0.1, col=2) nascente ### </code> ================================================ Éder Comunello Researcher at Brazilian Agricultural Research Corporation (Embrapa) DSc in Agricultural Systems Engineering (USP/Esalq) MSc in Environ. Sciences (UEM), Agronomist (UEM) --- Embrapa Agropecuária Oeste, Dourados, MS, Brazil |<O>| ================================================ GEO, -22.2752, -54.8182, 408m UTC-04:00 / DST: UTC-03:00 Em 9 de junho de 2016 08:41, Éder Comunello <comunello.eder@gmail.com> escreveu:
Cassiano, bom dia!
Imagino que o que você esteja tentando adicionar ao gráfico é uma linha da média + erro tal qual essa que aparece no gráfico abaixo:
Infelizmente, não encontrei como reproduzir esse efeito no R. O melhor que consegui foi marcar os limites com novas linhas. [image: Imagem inline 1]
### <code r> library(circular) nascente <- read.table(text= "AMBIENTE mes angle DISTANCIA ABUND RIQUEZA NASC 1 0 0 2 2 NASC 2 30 0 3 1 NASC 3 60 0 4 3 NASC 4 90 0 8 5 NASC 5 120 0 0 0 NASC 6 150 0 0 0 NASC 7 180 0 31 1 NASC 8 210 0 10 3 NASC 9 240 0 18 3 NASC 10 270 0 7 5 NASC 11 300 0 2 2 NASC 12 330 0 1 1", head=T)
nascente nascente0=nascente[which(nascente$DISTANCIA==0),]; nascente0 meses <- c("Mai","Jun","Jul","Ago","Set","Out","Nov","Dez","Jan","Fev","Mar","Abr")
#convertendo os dados para circular, usando radianos nasc.c0 <- circular(rep(rad(nascente0$angle),nascente0$ABUND), units=c("radians"), rota="clock") nasc.c0 plot(nasc.c0, stack=T)
# Angulo médio + Desvio padrão angular media <- mean(nasc.c0) dp <- sd(nasc.c0) n <- length(nasc.c0) err <- dp/sqrt(n) len <- rho.circular(nasc.c0) cbind(media, dp, n, err, len)
sd(nasc.c0); sqrt(-2*log(rho.circular(nasc.c0)))
# Gráfico # par(mfrow=c(1,3)) # par(mar=c(0,0,2,0)) plot(nasc.c0,axes=F,main=NA) axis.circular(at=circular(rad(nascente0$angle)), labels=meses, units=c("radians"), rota="clock") rose.diag(nasc.c0, add=T, axes=F, col=8, bor=8, lwd=3) # arrows.circular(median(nasc.c0),y=rho.circular(nasc.c0),lwd=3,lty=1,angle=30,len=0.1,col="green") #aqui eu inseri a mediana arrows.circular(mean(nasc.c0),y=len, lwd=2, lty=1, angle=90, len=0.1, col=2) arrows.circular(media+err, y=1, lwd=1, lty=3, angle=0, len=0.1, col=2) arrows.circular(media-err, y=1, lwd=1, lty=3, angle=0, len=0.1, col=2) ### </code>
================================================ Éder Comunello Researcher at Brazilian Agricultural Research Corporation (Embrapa) DSc in Agricultural Systems Engineering (USP/Esalq) MSc in Environ. Sciences (UEM), Agronomist (UEM) --- Embrapa Agropecuária Oeste, Dourados, MS, Brazil |<O>| ================================================ GEO, -22.2752, -54.8182, 408m UTC-04:00 / DST: UTC-03:00
participantes (2)
-
Cassiano -
Éder Comunello