# Carregando pacotes necessários:
library(tidyverse) # Impossível não precisar dos pacotes dentro do Tidyverse
library(gt) # Permite criar melhores tabelas e salvar em word ou pdf
library(gtExtras) # É uma extensão do pacote gt que fornece funções adicionais
library(gtsummary) # Utilizado para criar tabelas estatísticas de forma simples
library(labelled) # Funções para trabalhar com dados que contêm rótulos
library(gridExtra) # Organizar múltiplos gráficos em uma única página.
library(psych) # Faz análise descritiva robusta e sem muita programação
library(ggpubr) # Ajuda a encontrar valores da equação de regressão
library(patchwork) # Permite juntar gráficos na mesma janela
library(flextable) # Facilita o acesso para trablhar com tabelas
library(dlookr) # Ênfase em diagnóstico de qualidade dos dados em EDA
library(ggstatsplot)# Extensão do ggplot2 com gráficos de saídas + estatísticas
library(SmartEDA) # Simplificar e automatiza a EDA.
# Baixando o banco de dados em um dataframe direto do github
HRV <- read.csv(
"https://raw.githubusercontent.com/GleidsonUERN/LAFEX/main/VFC.csv")Análise Descritiva Exploratória
Realizando uma análise exploratória de dados (Exploratory Data Analysis - EDA).
Cenário
Nos dias atuais temos um vasto volume de informações que, se devidamente analisados, podem revelar padrões, tendências e insights valiosos que auxiliam na tomada de decisões estratégicas. No entanto, para extrair o máximo potencial desses dados, é fundamental realizar uma análise descritiva e exploratória minuciosa.
Análise descritiva.
A análise descritiva e exploratória de dados tem como principal objetivo resumir e interpretar as características principais de um conjunto de dados. Este tipo de análise é um passo essencial em qualquer projeto de ciência de dados, pois permite um entendimento profundo do cenário estudado, identificando comportamentos comuns, outliers, e relações entre variáveis.
Nesta análise, abordaremos um banco de dados específico (dados fictícios) mas que simulam dados originais. Exploreremos suas diversas através de técnicas estatísticas e ferramentas de visualização de dados, buscaremos responder perguntas como:
Quais são as principais características dos dados?
Como as variáveis se distribuem e quais são suas medidas de tendência central ou medidas de resumo (média, mediana, moda, etc.)?
Existem padrões ou tendências evidentes?
As variáveis dependentes são diferentes quando se toma uma variável categórica como independente?
Este processo não apenas facilita a compreensão do conjunto de dados, mas também orienta as etapas subsequentes de modelagem e predição. A análise descritiva e exploratória é, portanto, um pilar fundamental para qualquer estudo que vise transformar dados brutos em informações acionáveis e relevantes. Neste documento, apresentaremos os resultados dessa análise inicial, destacando os principais achados e preparando o terreno para análises mais aprofundadas. Convidamos você a embarcar conosco nesta jornada de descoberta e compreensão dos dados, onde cada gráfico, tabela e estatística, revelará um pouco mais sobre o universo que estamos investigando.
Carregando o dataset e os pacotes necessários
Certifique-se que você está online para fazer o download diretamente pelo script do R. Do contrário, faça o downloado do arquivo em https://github.com/GleidsonUERN/LAFEX/blob/main/VFC.csv e salve no seu diretório pessoal.
A seguir, temos os pacotes que poderão ser utilizados neste documento. Certifique-se que você tem todos eles instalados. Do contrário, use: install.packages("nome_pacote") para instalar àqueles que ainda não estão em sua máquina.
# Baixando o banco de dados em um dataframe direto do seu computador
HRV <- read.csv("C://pastas/subpastas/arquivo.csv")
Importante
Observe que a quantidade de pastas, subpastas bem como o disco rígido declarado na função read.csv() dependerá do local onde você salvou seu arquivo. Fique atento para passar a localização corretamente.
Conhecendo as variáveis do dataset
Agora vamos conhecer um pouco sobre o banco de dados e fazer algumas alterações.
# Exibindo as 6 primeiras linhas do dataframe
head(HRV) sex phys_Act age mRR SDNN rMSSD LF LFnu HF HFnu LF_HF
1 M Sedentário 28 836.38 75.22 47.92 632.49 37.52 1352.57 27.16 1.26
2 M Não Sedentário 27 832.40 55.15 28.46 193.53 45.03 82.82 19.74 2.34
3 M Não Sedentário 28 924.38 56.56 72.48 549.74 77.98 891.48 62.05 1.26
4 F Não Sedentário 30 914.10 49.83 46.35 342.29 51.63 1329.44 42.31 1.26
5 M Não Sedentário 37 785.21 65.38 26.62 991.16 61.13 2368.29 43.76 1.26
6 F Não Sedentário 30 1019.65 56.90 19.72 1079.52 50.15 362.07 28.07 2.98# Exibindo os nomes e os tipos das variáveis da base de dados.
dplyr::glimpse(HRV)Rows: 500
Columns: 11
$ sex <chr> "M", "M", "M", "F", "M", "F", "F", "F", "M", "M", "F", "F", "…
$ phys_Act <chr> "Sedentário", "Não Sedentário", "Não Sedentário", "Não Sedent…
$ age <int> 28, 27, 28, 30, 37, 30, 35, 33, 29, 23, 28, 28, 30, 36, 40, 3…
$ mRR <dbl> 836.38, 832.40, 924.38, 914.10, 785.21, 1019.65, 948.48, 1143…
$ SDNN <dbl> 75.22, 55.15, 56.56, 49.83, 65.38, 56.90, 32.83, 32.83, 46.95…
$ rMSSD <dbl> 47.92, 28.46, 72.48, 46.35, 26.62, 19.72, 25.91, 93.15, 27.98…
$ LF <dbl> 632.49, 193.53, 549.74, 342.29, 991.16, 1079.52, 1114.68, 377…
$ LFnu <dbl> 37.52, 45.03, 77.98, 51.63, 61.13, 50.15, 58.10, 51.47, 65.64…
$ HF <dbl> 1352.57, 82.82, 891.48, 1329.44, 2368.29, 362.07, 2482.94, 18…
$ HFnu <dbl> 27.16, 19.74, 62.05, 42.31, 43.76, 28.07, 51.33, 58.39, 55.28…
$ LF_HF <dbl> 1.26, 2.34, 1.26, 1.26, 1.26, 2.98, 1.26, 1.26, 1.26, 12.98, …
Repare que as variáveis sex e phys_Act estão alocadas como character. Desse modos, é importante que mudemos para uma variável categórica (factor). Existem inúmeras formas de se fazer isso. Vejamos abaixo uma forma bem elementar.
# Primeiro a variável sex
HRV$sex <- factor(HRV$sex)
# Agora vamos para a variável phys_Act
HRV$phys_Act <- factor(HRV$phys_Act)
# Verificando o resultado
dplyr::glimpse(HRV)Rows: 500
Columns: 11
$ sex <fct> M, M, M, F, M, F, F, F, M, M, F, F, F, M, F, M, F, M, M, M, M…
$ phys_Act <fct> Sedentário, Não Sedentário, Não Sedentário, Não Sedentário, N…
$ age <int> 28, 27, 28, 30, 37, 30, 35, 33, 29, 23, 28, 28, 30, 36, 40, 3…
$ mRR <dbl> 836.38, 832.40, 924.38, 914.10, 785.21, 1019.65, 948.48, 1143…
$ SDNN <dbl> 75.22, 55.15, 56.56, 49.83, 65.38, 56.90, 32.83, 32.83, 46.95…
$ rMSSD <dbl> 47.92, 28.46, 72.48, 46.35, 26.62, 19.72, 25.91, 93.15, 27.98…
$ LF <dbl> 632.49, 193.53, 549.74, 342.29, 991.16, 1079.52, 1114.68, 377…
$ LFnu <dbl> 37.52, 45.03, 77.98, 51.63, 61.13, 50.15, 58.10, 51.47, 65.64…
$ HF <dbl> 1352.57, 82.82, 891.48, 1329.44, 2368.29, 362.07, 2482.94, 18…
$ HFnu <dbl> 27.16, 19.74, 62.05, 42.31, 43.76, 28.07, 51.33, 58.39, 55.28…
$ LF_HF <dbl> 1.26, 2.34, 1.26, 1.26, 1.26, 2.98, 1.26, 1.26, 1.26, 12.98, …
Note que as 2 variáveis agora são do tipo fator (factor), cada uma com 2 levels (categorias). Isso será muito importante para definirmos os grupos que irão confrontar as variáveis contínuas.
Plotando os primeiros gráficos
(Inspeção visual da distribuição)
Agora vamos plotar alguns gráficos do tipo histograma para visualizarmos a distribuição de algumas das variáveis contínuas.
# Primeiro a variável "age" (Idade)
plot_age <- hist(HRV$age,
col = "blue",
xlab = "Idade em anos",
main = "Histograma",
freq = F)
# Agora com a variável mRR (média dos intervalos R-R)
hist(HRV$mRR,
col = "red",
xlab = "Média dos Intervalos R-R em ms",
main = "Histograma",
freq = F)
Usando a função hist, nativa do R, podemos ir gerando um a um, todos os histogramas de cada uma das variáveis. Observe que não foi separado por grupo (Ex.: sex) e todos os casos da variável estão juntos. Teríamos que adicionar também a função curve(dnorm(x, mean, sd, add)) em adição à cada um dos gráficos. Isso não é nada amigável.
Podemos usar o pacote ggplot2 para criar algo mais eficiente e elegante. Vamos criar um histograma usando a função geom_histogram(), vamos adicionar a curva normal usando a função stat_function(), depois vamos criar um boxplot em cima do histograma usando a função geom_boxplot(). A função stat_boxplot() apenas adicionará linhas transversais nas extremidades dos bigodes do boxplot (barra de erros) e a função facet_grid()irá separar as plotagens (figuras) em função do sexo.
library(ggplot2)
ggplot(HRV, aes(x = age)) +
geom_histogram(aes(x = age, y = after_stat(density)),
fill = "lightblue",
color = "black",
bins = max(HRV$age) - min(HRV$age)) +
stat_function(fun = dnorm, args = list(mean = mean(HRV$age), sd = sd(HRV$age))) +
geom_boxplot(data = HRV, aes(x = age),
color = "red",
size = 0.9,
width = max(density(HRV$age)$y)/3,
outlier.colour = "red",
outlier.size = 4) +
stat_boxplot(geom = "errorbar", width = max(density(HRV$age)$y)/5) +
facet_grid(sex ~ .)
Veja como ficou elegante. Temos um histograma plotado abaixo de um boxplot da mesma variável e do mesmo grupo. Com esse tipo de visualização já é possível ver que temos um outlier na idade do grupo Masculino. Ainda assim, visualmente podemos crer que a distribuição destas 2 variáveis respeitam os limites da curva normal para que possamos chamá-las de paramétrica.
Otimizando a saída dos gráficos para datasets grandes
Pense que independente de usarmos a função hist() nativa do R ou funções mais densas como as do ggplot(), em princípio teríamos que passar todas as variáveis manualmente para criarmos todos os gráficos. Então, porque não otimizarmos isso para que com um pouco mais de esforço (entendimento do código), tenhamos um resultado melhor.
Aviso
Se você não está muito disposto de avançar na escrita de código e resolução de problemas e deseja continuar plotando os gráficos 1 à 1, então é melhor pular este tópico.
Faremos agora uma função que irá plotar todos os gráficos que precisamos do nosso dataframe. O que faremos é usar o trecho de código acima como parte principal da função. Adicionaremos um argumento que levará em conta todos os nomes das variáveis numéricas do dataframe. Também irei modificar o argumento aes() para aes_string() já que utilizaremos uma string para declarar o nome da variável. Tente ler o código e pensar o que ele está executanto. Se tiver dificuldades de reproduzir este código, aplicando às características do seu dataframe, use o trecho de código anterior modificando manualmente o nome das variáveis (uma por vez) de acordo com o seu dataframe.
Bem, vamos criar a função que executará repetidamente a criação do gráfico mostrado anteriormente baseado em 1 parâmetro que precisaremos estabelecer. O parâmetro var_name usará os nomes das variáveis do dataframe para rodar a função adequadamente. Aqui neste exemplo, continuarei usando o nosso dataframe HRV, então note que várias passagens do código já está inserido o nome do dataframe HRV e da variável sex. Desse modo, ao tentar copiar o código e utilizar em seu próprio dataframe, você precisará modificar no código da função: a fonte dos dados e a variável que irá realizar a divisão pelo facet_grid().
# Unindo várias funções ggplot para criar um gráfico que represente vários
# aspectos da distribuição da variável.
criando_graficos <- function(var_name) {
# !!sym foi necessário pois o argumento precisará ser dinâmico
ggplot(HRV, aes(x = !!sym(var_name))) +
# Aqui temos o histograma
geom_histogram(aes(y = after_stat(density)),
fill = "lightblue",
color = "black",
bins = 30) +
# Aqui temos a adição da curva normal
stat_function(fun = dnorm,
args = list(mean = mean(HRV[[var_name]],
na.rm = TRUE),
sd = sd(HRV[[var_name]],
na.rm = TRUE)),
color = "blue", linewidth = 1) +
# Aqui temos a adição do boxplot com parâmetros que se adequará à cada
# amplitude das variáveis.
geom_boxplot(aes(y = 0),
color = "red",
size = 0.9,
width = max(density(HRV[[var_name]])$y)/3,
outlier.colour = "red",
outlier.size = 3) +
# Aqui estamos adicionando a barra de erros no boxplot
stat_boxplot(geom = "errorbar", width = max(density(HRV[[var_name]])$y)/5) +
# Aqui estamos fazendo a plotagem em função da variável "sex"
facet_grid(sex ~ .) +
labs(x = var_name, y = "Densidade", title = paste("Distribuição de", var_name))
}É importante lembrar que a função acima foi criada por nós e, portanto, não se refere à nenhum pacote específico. Ela permanecerá no ambiente global enquanto o script estiver ativo. Considerando nosso dataframe HRV, já sabemos que ele tem variáveis que são numéricas e categóricas. Para a criação de histogramas e boxplots, somente as variáveis numéricas nos interessam. Vamos usar a função sapply() para guardar apenas as variáveis numéricas.
# Selecionar as variáveis numéricas
numeric_vars <- names(HRV)[sapply(HRV, is.numeric)]
print(numeric_vars)[1] "age" "mRR" "SDNN" "rMSSD" "LF" "LFnu" "HF" "HFnu" "LF_HF"
Note que a variável criada (numeric_vars) abriga os nomes das variáveis do dataframe HRV que são numéricas. Os nomes estão em formato de string (por isso utilizamos !!sym() no argumento do ggplot na função que criamos).
Um pouco sobre o pacote
purrr antes de continuar
Imagine agora que você tem uma lista de coisas para fazer e quer aplicar a mesma ação em cada item dessa lista. Ou seja, temos uma lista de variáveis e queremos aplicar a função de ciação do gráfico em cada uma dessas variáveis. Fazer isso manualmente seria bem chato, né?
Aí que entra a função map() do pacote purrr. Ela é tipo um super assistente que pega essa sua lista de coisas e aplica uma função em cada item dessa lista, de um jeito muito mais rápido e fácil. Ele vai lá, pega cada número/string da lista, posiciona no local determinado dentro da função e te devolve uma nova lista com tudo já bonitinho.
O legal é que você pode usar o map() para fazer várias coisas, não só com números, mas com listas de textos, data frames, ou o que você imaginar. Quer transformar uma lista de nomes para maiúsculas? O map() faz isso pra você. Quer aplicar uma função que você criou em cada coluna de um data frame? O map() resolve. Resumindo: o map() é tipo aquele amigo super eficiente que te ajuda a aplicar uma função em cada item de uma lista de forma rápida e sem complicação. Então vamos lá!
Como criamos nossa função anteriormente, então agora é hora de aplicar a nossa função como um argumento dentro da função map() juntamente com o vetor (numeric_vars) que será passado linearmente dentro da função.
# Criando os gráficos para todas as variáveis numéricas
meus_plots <- purrr::map(numeric_vars, criando_graficos)
# A função map() retorna uma lista. No nosso caso foi uma lista com 9 plots
# Então vamos exibí-los em uma única figura com o seguinte layout (3 x 3). Fique
# atento a quantidade de variáveis que voçê tem para que possa modificar
# corretamente a quantidade de linhas (nrow) e colunas (ncol).
grid.arrange(grobs = meus_plots, nrow = 3, ncol = 3)
Outras alternativas com os pacotes gt e gtExtras.
Se você chegou até aqui, deve estar apresentando algum nível de cansaço ou está avaliando se o resultado gráfico exibido até agora justificam tantas linhas de código. Bem, a ideia é fazer você perceber a versatilidade e plasticidade que temos para criar gráficos da maneira que bem queremos. Algo bem diferente dos gráficos clássicos feitos em Excel ou mesmo em GraphPad. Contudo, se você é do tipo que está pouco se importando com algum requinte gráfico para exibição dos dados então aí vai um alento: veja o que a função gt_plt_summary() do pacote gtExtras pode fazer.
HRV |> gtExtras::gt_plt_summary(title = "Análise descritiva da base de dados HRV")
Note que a saída nos informou diversas características da base de dados (dataframe). Na parte esquerda, temos ícones que descrevem o tipo da variável. Podemos observar que há uma diferença entre as variáveis categóricas e as variáveis contínuas. Em seguida o nome da variável conforme consta no data frame. Depois temos um miniplot em forma de histograma para cada uma das variáveis com destaque para as variáveis contínuas que apresentam um pequeno histograma com a média posicionada (linha vertical de cor preta) em cada uma delas. Podemos ainda verificar dados ausentes (missings) e algumas medidas tais como média (MEAN), mediana (MEDIAN) e desvio padrão (SD). Para quem trabalha com base de dados grandes, este comando é muito eficiente para que tenhamos uma visão geral das variáveis que vamos trabalhar sem muito esforço.
Visualizando informações em tabelas
Vamos agora considerar uma visualização em tabela para exploração dos dados. Muitos são os pacotes envolvidos nestas tarefas. Pacotes como o psych que tem uma função describe() (semelhante à mesma função do pacote dlookr). Cada um, com suas potencialidades e limitações. Aqui em especial utilizaremos a função describe() do pacote dlookr.
# O pacote `dlookr` desconsidera as variáveis categóricas na análise.
# Isso pode ser interessante quando não tempos grupos a serem considerados.
dlookr::describe(HRV) |>
dplyr::select(described_variables, mean, sd, se_mean, skewness, kurtosis) |>
flextable::flextable()described_variables | mean | sd | se_mean | skewness | kurtosis |
|---|---|---|---|---|---|
age | 30.11000 | 4.539064 | 0.2029931 | -0.13430470 | -0.051492409 |
mRR | 928.73874 | 88.535808 | 3.9594417 | 0.08666669 | 0.007503591 |
SDNN | 50.63884 | 14.219739 | 0.6359261 | 0.50400889 | -0.461228198 |
rMSSD | 42.80574 | 19.685768 | 0.8803743 | 0.72543429 | 0.024277849 |
LF | 539.99948 | 246.192526 | 11.0100645 | 0.47005279 | -0.398464768 |
LFnu | 51.86088 | 9.727337 | 0.4350197 | 0.04795421 | -0.169713146 |
HF | 735.64148 | 613.161768 | 27.4214279 | 0.76203099 | -0.227598678 |
HFnu | 39.54764 | 9.981022 | 0.4463649 | 0.17960129 | 0.246069653 |
LF_HF | 2.69504 | 2.742045 | 0.1226280 | 2.08759979 | 3.763881690 |
# Vamos salvar as informações em um dataframe
minha_EDA <- dlookr::describe(HRV) |>
dplyr::select(described_variables, mean, sd, se_mean, skewness, kurtosis) |>
dplyr::mutate("Coef. Variação" = (sd/mean)*100)
A saída da função describe()do pacote dlookr tem muitas colunas (cada uma delas com um tipo de dado estatístico diferente). Desse modo, utilizamos a função select() do dplyr para escolher apenas estatísticas de interesse. Foram elas: as próprias variaéis (described_variables), a média (mean), desvio padrão (sd) a assimetria (skewness) e o achatamento (kurtosis). Se você passar a função sozinha (dlookr::describe(HRV)) verá o que isso quer dizer.
Conhecendo um pouco sobre Assimetria e Achatamento
Skewness, ou assimetria, é uma medida que descreve a simetria dos dados em torno da média. Ela indica se a distribuição dos dados é simétrica ou inclinada para algum dos lados. Existem três tipos principais de skewness como veremos a seguir.
Assimetria Positiva (Right-Skewed)
- A cauda direita da distribuição é mais longa ou “gorda” do que a cauda esquerda. A maioria dos valores está concentrada no lado esquerdo, com alguns valores extremos à direita. Exemplo: Renda de uma população onde poucas pessoas ganham muito mais do que a maioria.
Assimetria Negativa (Left-Skewed)
- A cauda esquerda da distribuição é mais longa ou “gorda” do que a cauda direita. A maioria dos valores está concentrada no lado direito, com alguns valores extremos à esquerda. Exemplo: Idade de aposentadoria onde poucas pessoas se aposentam muito antes do comum.
Assimetria Nula (Symmetric)
- A distribuição é simétrica em torno da média. A cauda direita e a cauda esquerda têm aproximadamente o mesmo comprimento. Exemplo: Altura de uma população onde os valores se distribuem igualmente ao redor da média.
A skewness é importante na análise de dados porque uma distribuição assimétrica pode afetar muitas técnicas estatísticas que assumem uma distribuição normal. Compreender a skewness ajuda a decidir se é necessário aplicar transformações nos dados para aproximá-los de uma distribuição de referência como a normal.
Nota
Assimetria (Skewness)
Existem limites estabelecidos por diferentes autores que ajudam a identificar quando uma distribuição pode ser considerada normal. Aqui estão os valores geralmente aceitos:
Valores de Assimetria (Skewness)
Assimetria entre -0.5 e 0.5: Considerada simétrica.
Assimetria entre -1 e -0.5 ou entre 0.5 e 1: Levemente assimétrica.
Assimetria menor que -1 ou maior que 1: Fortemente assimétrica.
Curtose, por sua vez, é uma medida que descreve a forma das caudas de uma distribuição em comparação com uma distribuição normal. Ela indica a presença de outliers e a “gordura” das caudas. Existem três tipos principais de kurtosis:
Leptocúrtica (Leptokurtic)
- Distribuições com caudas mais pesadas ou “gordas” e um pico mais alto e estreito do que a distribuição normal. Isso indica a presença de muitos outliers. Exemplo: Notas de um exame onde a maioria dos alunos tira notas muito altas ou muito baixas, mas poucos ficam na média.
Mesocúrtica (Mesokurtic)
- Distribuições com caudas e pico semelhantes aos da distribuição normal. Exemplo: Altura de uma população, se a distribuição for semelhante à distribuição normal.
Platicúrtica (Platykurtic)
- Distribuições com caudas mais leves ou “finas” e um pico mais baixo e largo do que a distribuição normal. Isso indica menos outliers. Exemplo: Peso de uma população onde a maioria dos valores está concentrada em torno da média, com poucas variações extremas.
Nota
Kurtosis (Curtose ou Achatamento)
A kurtosis é importante na análise de dados porque indica a probabilidade de ocorrência de valores extremos (outliers). Distribuições com alta kurtosis (leptocúrticas) têm maior probabilidade de ter outliers, o que pode afetar a robustez de várias análises estatísticas. Aqui estão os valores geralmente aceitos:
Valores de Curtose (Kurtosis)
A curtose excessiva é medida comparando-se a curtose da distribuição dos dados com a curtose da distribuição normal (que é 3).
Curtose entre -1 e 1: Considerada aceitável para a normalidade.
Curtose entre -2 e 2: Muitas vezes ainda considerada razoável.
Curtose menor que -2 ou maior que 2: Indicativa de uma distribuição não normal, com caudas muito leves ou muito pesadas (Curtose excessiva).
Os limites para considerar uma distribuição normal variam ligeiramente entre diferentes autores, mas a faixa de assimetria entre -1 e +1 e de curtose excessiva entre -2 e +2 é amplamente aceita na literatura como razoável para indicar normalidade. Esses parâmetros ajudam a orientar a interpretação de dados e a escolha de transformações quando necessário (Bulmer 1979; Hair 2009; Kline 2023; Mallery e George 2000; Ptaszyński 2022).
Baseado nisso podemos criar uma coluna para classificar os valores de Assimetria e Curtose que temos até agora e tomar alguma decisão.
# Criando uma função para classificar as variáveis com base em skewness e kurtosis
classificacao_variavel <- function(skewness, kurtosis) {
if (skewness >= -1 & skewness <= 1 & kurtosis >= -2 & kurtosis <= 2) {
return("Paramétrica")
} else {
return("Não-Paramétrica")
}
}
# Aplicando a função para criar a coluna de classificação
minha_EDA <- minha_EDA %>%
mutate(classificacao = mapply(classificacao_variavel, skewness, kurtosis))
minha_EDA |> flextable::flextable()described_variables | mean | sd | se_mean | skewness | kurtosis | Coef. Variação | classificacao |
|---|---|---|---|---|---|---|---|
age | 30.11000 | 4.539064 | 0.2029931 | -0.13430470 | -0.051492409 | 15.074940 | Paramétrica |
mRR | 928.73874 | 88.535808 | 3.9594417 | 0.08666669 | 0.007503591 | 9.532908 | Paramétrica |
SDNN | 50.63884 | 14.219739 | 0.6359261 | 0.50400889 | -0.461228198 | 28.080697 | Paramétrica |
rMSSD | 42.80574 | 19.685768 | 0.8803743 | 0.72543429 | 0.024277849 | 45.988618 | Paramétrica |
LF | 539.99948 | 246.192526 | 11.0100645 | 0.47005279 | -0.398464768 | 45.591252 | Paramétrica |
LFnu | 51.86088 | 9.727337 | 0.4350197 | 0.04795421 | -0.169713146 | 18.756599 | Paramétrica |
HF | 735.64148 | 613.161768 | 27.4214279 | 0.76203099 | -0.227598678 | 83.350624 | Paramétrica |
HFnu | 39.54764 | 9.981022 | 0.4463649 | 0.17960129 | 0.246069653 | 25.237971 | Paramétrica |
LF_HF | 2.69504 | 2.742045 | 0.1226280 | 2.08759979 | 3.763881690 | 101.744137 | Não-Paramétrica |
De acordo com as condições que foram estabelecidas, apenas a variável LF_HF se apresentou como uma variável não-paramétrica. Isso implicará em uma utilização de teste estatístico diferente do teste utilizado com as variáveis paramétricas. Note ainda que a função describe() tomou como conjunto de dados, todas as entradas do dataframe e não as separou por nenhuma variável categórica. Desse modo, temos que calcular as mesmas medidas levando em conta cada variável categórica que utilizaremos em uma comparação. Atualizarei, no futuro, este script com funções que permitirão exibir a tabela separado por determinadas variáveis de interesse (variáveis categóricas). Outras funções do pacote dlookr merecem ser mencionadas para o uso frequente desse pacote.
# Verificar o tipo das variáveis e dados ausentes
dlookr::diagnose(HRV) |> flextable::flextable()variables | types | missing_count | missing_percent | unique_count | unique_rate |
|---|---|---|---|---|---|
sex | factor | 0 | 0 | 2 | 0.004 |
phys_Act | factor | 0 | 0 | 2 | 0.004 |
age | integer | 0 | 0 | 26 | 0.052 |
mRR | numeric | 0 | 0 | 497 | 0.994 |
SDNN | numeric | 0 | 0 | 416 | 0.832 |
rMSSD | numeric | 0 | 0 | 397 | 0.794 |
LF | numeric | 0 | 0 | 445 | 0.890 |
LFnu | numeric | 0 | 0 | 459 | 0.918 |
HF | numeric | 0 | 0 | 381 | 0.762 |
HFnu | numeric | 0 | 0 | 461 | 0.922 |
LF_HF | numeric | 0 | 0 | 192 | 0.384 |
# Análise descritiva mantendo as variáveis contínuas originais em linha
dlookr::diagnose_numeric(HRV) |> flextable::flextable()variables | min | Q1 | mean | median | Q3 | max | zero | minus | outlier |
|---|---|---|---|---|---|---|---|---|---|
age | 17.00 | 27.0000 | 30.11000 | 30.000 | 33.0000 | 42.00 | 0 | 0 | 1 |
mRR | 712.03 | 864.5575 | 928.73874 | 931.440 | 985.4350 | 1,231.13 | 0 | 0 | 4 |
SDNN | 30.35 | 38.0750 | 50.63884 | 50.025 | 59.6225 | 92.31 | 0 | 0 | 1 |
rMSSD | 18.76 | 25.2125 | 42.80574 | 40.300 | 56.0175 | 117.55 | 0 | 0 | 2 |
LF | 191.28 | 350.4300 | 539.99948 | 512.880 | 713.4925 | 1,319.57 | 0 | 0 | 1 |
LFnu | 29.63 | 45.8175 | 51.86088 | 51.405 | 58.2025 | 78.54 | 0 | 0 | 3 |
HF | 80.08 | 105.0300 | 735.64148 | 621.660 | 1,148.9875 | 2,658.77 | 0 | 0 | 0 |
HFnu | 15.11 | 32.8175 | 39.54764 | 39.635 | 46.0700 | 74.46 | 0 | 0 | 4 |
LF_HF | 0.92 | 1.2600 | 2.69504 | 1.260 | 2.4700 | 14.77 | 0 | 0 | 94 |
# Análise descritiva mantendo as variáveis categóricas originais em linha
dlookr::diagnose_category(HRV) |> flextable::flextable()variables | levels | N | freq | ratio | rank |
|---|---|---|---|---|---|
sex | M | 500 | 261 | 52.2 | 1 |
sex | F | 500 | 239 | 47.8 | 2 |
phys_Act | Sedentário | 500 | 254 | 50.8 | 1 |
phys_Act | Não Sedentário | 500 | 246 | 49.2 | 2 |
# Verificando os outliers e gerando uma tabela flextable
dlookr::diagnose_outlier(HRV) |> flextable::flextable()variables | outliers_cnt | outliers_ratio | outliers_mean | with_mean | without_mean |
|---|---|---|---|---|---|
age | 1 | 0.2 | 17.000000 | 30.11000 | 30.136273 |
mRR | 4 | 0.8 | 1,200.590000 | 928.73874 | 926.546391 |
SDNN | 1 | 0.2 | 92.310000 | 50.63884 | 50.555331 |
rMSSD | 2 | 0.4 | 113.965000 | 42.80574 | 42.519960 |
LF | 1 | 0.2 | 1,319.570000 | 539.99948 | 538.437214 |
LFnu | 3 | 0.6 | 78.053333 | 51.86088 | 51.702777 |
HF | 0 | 0.0 | 735.64148 | 735.641480 | |
HFnu | 4 | 0.8 | 71.535000 | 39.54764 | 39.289677 |
LF_HF | 94 | 18.8 | 7.784043 | 2.69504 | 1.516798 |
Destaco a função diagnose_outlier() que nos apresenta o total de outliers, a frequência que eles representam do conjunto de dados e o cálculo da média COM (with) e SEM (without) os outliers. Podemos até plotar um gráfico para visualizar essas informações. Vejamos.
library(dlookr)
plot_outlier(HRV)
Outras alternativas para visualização em tabelas incluem a família dos pacotes relacionados ao pacote gt tais como o gtExtras e o gtsummary. Estes 3 pacotes permitem uma infinidade de funções para se trabalhar com tabelas descritivas e inferenciais.
HRV |>
gtsummary::tbl_summary(by = sex) |>
gtsummary::add_overall(last = T)| Characteristic | F, N = 2391 | M, N = 2611 | Overall, N = 5001 |
|---|---|---|---|
| phys_Act | |||
| Não Sedentário | 118 (49%) | 128 (49%) | 246 (49%) |
| Sedentário | 121 (51%) | 133 (51%) | 254 (51%) |
| age | 31 (28, 34) | 30 (27, 33) | 30 (27, 33) |
| mRR | 936 (863, 984) | 930 (867, 986) | 931 (865, 985) |
| SDNN | 50 (39, 60) | 50 (37, 59) | 50 (38, 60) |
| rMSSD | 40 (23, 56) | 40 (27, 56) | 40 (25, 56) |
| LF | 519 (378, 711) | 506 (336, 712) | 513 (350, 713) |
| LFnu | 51 (45, 58) | 52 (46, 59) | 51 (46, 58) |
| HF | 554 (116, 1,176) | 640 (83, 1,114) | 622 (105, 1,149) |
| HFnu | 40 (33, 46) | 39 (33, 46) | 40 (33, 46) |
| LF_HF | 1.26 (1.26, 2.64) | 1.26 (1.26, 2.34) | 1.26 (1.26, 2.47) |
| 1 n (%); Median (IQR) | |||
Com o pacote labelled podemos criar labels para os nomes das variáveis do dataframe. Desse modo, a visualização de informações em gráficos e tabelas utilizará os labels em lugar do nome original. Veja abaixo.
# Vamos salvar os labels no dataframe
HRV <- HRV %>%
labelled::set_variable_labels(sex = "Sexo",
phys_Act = "Nível de Atividade Física",
age = "Idade (anos)",
LF_HF = "LF/HF")
# Reseta definições de tema
reset_gtsummary_theme()
# Traduz informações estatísticas para o Port.
theme_gtsummary_language("pt")
HRV |>
gtsummary::tbl_summary(by = sex) |>
gtsummary::add_overall(last = T)| Características | F, N = 2391 | M, N = 2611 | Total, N = 5001 |
|---|---|---|---|
| Nível de Atividade Física | |||
| Não Sedentário | 118 (49%) | 128 (49%) | 246 (49%) |
| Sedentário | 121 (51%) | 133 (51%) | 254 (51%) |
| Idade (anos) | 31 (28, 34) | 30 (27, 33) | 30 (27, 33) |
| mRR | 936 (863, 984) | 930 (867, 986) | 931 (865, 985) |
| SDNN | 50 (39, 60) | 50 (37, 59) | 50 (38, 60) |
| rMSSD | 40 (23, 56) | 40 (27, 56) | 40 (25, 56) |
| LF | 519 (378, 711) | 506 (336, 712) | 513 (350, 713) |
| LFnu | 51 (45, 58) | 52 (46, 59) | 51 (46, 58) |
| HF | 554 (116, 1,176) | 640 (83, 1,114) | 622 (105, 1,149) |
| HFnu | 40 (33, 46) | 39 (33, 46) | 40 (33, 46) |
| LF/HF | 1.26 (1.26, 2.64) | 1.26 (1.26, 2.34) | 1.26 (1.26, 2.47) |
| 1 n (%); Mediana (AIQ) | |||
O default dessa função utiliza a mediana. Vamos usar um tema que utiliza o desvio padrão e a média.
gtsummary::reset_gtsummary_theme()
gtsummary::theme_gtsummary_language("pt")
gtsummary::theme_gtsummary_mean_sd()
HRV |>
gtsummary::tbl_summary(by = sex, digits = everything() ~ 2) |>
gtsummary::add_overall(last = T)| Características | F, N = 2391 | M, N = 2611 | Total, N = 5001 |
|---|---|---|---|
| Nível de Atividade Física | |||
| Não Sedentário | 118.00 (49.37%) | 128.00 (49.04%) | 246.00 (49.20%) |
| Sedentário | 121.00 (50.63%) | 133.00 (50.96%) | 254.00 (50.80%) |
| Idade (anos) | 30.40 (4.47) | 29.85 (4.60) | 30.11 (4.54) |
| mRR | 928.63 (93.35) | 928.84 (84.06) | 928.74 (88.54) |
| SDNN | 51.07 (14.45) | 50.24 (14.03) | 50.64 (14.22) |
| rMSSD | 42.63 (20.03) | 42.97 (19.40) | 42.81 (19.69) |
| LF | 546.19 (241.24) | 534.34 (250.97) | 540.00 (246.19) |
| LFnu | 51.39 (9.74) | 52.29 (9.72) | 51.86 (9.73) |
| HF | 732.92 (624.89) | 738.13 (603.41) | 735.64 (613.16) |
| HFnu | 39.67 (10.06) | 39.44 (9.93) | 39.55 (9.98) |
| LF/HF | 2.73 (2.76) | 2.66 (2.73) | 2.70 (2.74) |
| 1 n (%); Média (Desvio Padrão) | |||
Testes específicos para análise da distribuição dos dados.
Até agora, foi possível conhecer vários aspectos das variáveis de nosso dataset. Tanto pela inspeção visual dos gráficos como pelas medidas de ajuste da curva de distribuição dos dados é possível tomar algumas decisões. Apesar disso, nos dias atuais, muitas revistas e revisores querem se deparar com testes que apoiem a classificação da distribuição dos seus dados. Então é hora de ir além! Vamos investigar objetivamente a distribuição das variáveis com testes de hipóteses para tal fim. Os testes de Shapiro-Wilk e Kolmogorov-Smirnov são ambos usados para avaliar a normalidade de uma distribuição de dados, mas possuem diferenças significativas em termos de aplicação e funcionamento (Shapiro e Wilk 1965; Lilliefors 1967; Moore, McCabe, e Craig 2009).
Nota
Teste de Shapiro-Wilk
O teste calcula uma estatística que compara a amostra de dados com uma distribuição normal ideal. Ele mede a correlação entre os dados da amostra e os valores normais teóricos esperados. Quanto maior a correlação, mais próxima a distribuição da amostra está de uma distribuição normal.
É mais adequado para amostras menores (há muita divergência sobre isso mas, a maioria dos autores concordam que a sensibilidade do teste é boa com até 50 amostras) (Gupta et al. 2019).
Interpretação:
- Hipótese nula (H0): A amostra segue uma distribuição normal.
- Hipótese nula (H1): Se o valor p resultante do teste for menor que um nível de significância escolhido (por exemplo, 0,05), rejeitamos a hipótese nula, indicando que a amostra não segue uma distribuição normal.
Nota
Teste de Kolmogorov-Smirnov (K-S)
Compara uma amostra de dados com uma distribuição de referência (por exemplo, normal, exponencial) ou comparar duas amostras entre si. Baseia-se na distância máxima (estatística D) entre a função de distribuição acumulada (CDF) empírica da amostra e a CDF da distribuição de referência. Essa distância é usada para avaliar se a amostra difere significativamente da distribuição de referência (frequentemente se usa a distribuição normal como sendo referência). Pode ser usado para qualquer tipo de distribuição de referência, não se limitando à normalidade.
Adequado para amostras maiores que 50 unidades e também para comparar duas amostras. Menos sensível em detectar desvios de normalidade em amostras pequenas comparado ao teste de Shapiro-Wilk (Gupta et al. 2019).
Interpretação
- Hipótese nula (H0): A distribuição da amostra corresponde à distribuição de referência (ou as duas amostras vêm da mesma distribuição).
- Hipótese alternativa (H1): Se o valor p resultante do teste for menor que um nível de significância escolhido (por exemplo, 0,05), rejeitamos a hipótese nula, indicando que a distribuição da amostra difere da distribuição de referência (ou as duas amostras não vêm da mesma distribuição).
Como o tamanho amostral no dataset HRV é superior a 50 unidades em qualquer saída de variável dependente em função de uma variável categórica. Apesar da amostra total ser n = 500, temos que analisar as distribuições dos tamanhos que irão compor os grupos ao qual se formam de acordo com as variáveis independentes (geralmente categóricas). Desse modo, vamos utilizar a variável sexo (sex) como variável independente e analisar todas as outras variáveis exceto a variável phys_Act (Nível de Atividade Física) que também é uma variável categórica.
# Calculando o K-S teste no sexo Feminino. O resultado da função ks.teste será
# uma lista com os resultados do teste.
ks_fem <- HRV |> dplyr::filter(sex == "F") |>
dplyr::select(where(is.numeric)) |>
purrr::map(~ ks.test(.x, "pnorm", mean = mean(.x), sd = sd(.x)))
# Calculando o K-S teste no sexo Masculino
ks_masc <- HRV |> dplyr::filter(sex == "M") |>
dplyr::select(where(is.numeric)) |>
purrr::map(~ ks.test(.x, "pnorm", mean = mean(.x), sd = sd(.x)))
# Guardando as informações relevantes da lista em um dataframe
library(tibble)
df_fem <- tibble("Variáveis" = names(ks_fem),
"Estatística" = round(purrr::map_dbl(ks_fem, ~ .x$statistic), 4),
"P valor" = round(purrr::map_dbl(ks_fem, ~ .x$p.value), 4))
df_masc <- tibble("Variáveis" = names(ks_masc),
"Estatística" = round(purrr::map_dbl(ks_masc, ~ .x$statistic), 4),
"P valor" = round(purrr::map_dbl(ks_masc, ~ .x$p.value), 4))
# Criando um classificador do resultado com base em critérios estabelecidos
df_fem <- df_fem |>
dplyr::mutate(
"Classificação Fem" = dplyr::case_when(`P valor` <= 0.05 ~ "Não Paramétrico",
TRUE ~ "Paramétrico")) |>
print()# A tibble: 9 × 4
Variáveis Estatística `P valor` `Classificação Fem`
<chr> <dbl> <dbl> <chr>
1 age 0.0921 0.0345 Não Paramétrico
2 mRR 0.046 0.694 Paramétrico
3 SDNN 0.0777 0.112 Paramétrico
4 rMSSD 0.118 0.0026 Não Paramétrico
5 LF 0.0913 0.0373 Não Paramétrico
6 LFnu 0.0453 0.712 Paramétrico
7 HF 0.148 0.0001 Não Paramétrico
8 HFnu 0.0383 0.875 Paramétrico
9 LF_HF 0.284 0 Não Paramétrico df_masc <- df_masc |>
dplyr::mutate(
"Classificação Masc" = dplyr::case_when(`P valor` <= 0.05 ~ "Não Paramétrico",
TRUE ~ "Paramétrico")) |>
print()# A tibble: 9 × 4
Variáveis Estatística `P valor` `Classificação Masc`
<chr> <dbl> <dbl> <chr>
1 age 0.0566 0.372 Paramétrico
2 mRR 0.036 0.888 Paramétrico
3 SDNN 0.087 0.0385 Não Paramétrico
4 rMSSD 0.109 0.0039 Não Paramétrico
5 LF 0.0861 0.0416 Não Paramétrico
6 LFnu 0.0394 0.813 Paramétrico
7 HF 0.139 0.0001 Não Paramétrico
8 HFnu 0.0298 0.974 Paramétrico
9 LF_HF 0.312 0 Não Paramétrico # Exibindo um resumo com as informações do teste
resumo <- tibble::tibble("Variáveis Dependentes" = df_fem$Variáveis,
"Amostra Feminina" = df_fem$`Classificação Fem`,
"Amostra Masculina" = df_masc$`Classificação Masc`)
resumo |> gt::gt() |>
gt::cols_align(align = c("center"))| Variáveis Dependentes | Amostra Feminina | Amostra Masculina |
|---|---|---|
| age | Não Paramétrico | Paramétrico |
| mRR | Paramétrico | Paramétrico |
| SDNN | Paramétrico | Não Paramétrico |
| rMSSD | Não Paramétrico | Não Paramétrico |
| LF | Não Paramétrico | Não Paramétrico |
| LFnu | Paramétrico | Paramétrico |
| HF | Não Paramétrico | Não Paramétrico |
| HFnu | Paramétrico | Paramétrico |
| LF_HF | Não Paramétrico | Não Paramétrico |
Observe que temos uma variedade de situações baseadas nas classficações das distribuições. Algumas das comparações poderão ser feitas considerando que ambos os grupos apresentaram distribuição normal.Desse modo, o teste estatístico poderá ser da família t (escolhendo adequadamente de acordo com o tipo da comparação e se as variâncias forem iguais). Contudo, se um ou os dois grupos apresenta uma distribuição não-paramétrica, então teremos que utilizar um teste equivalente não-paramétrico para realizar o teste de comparação.
Escolhendo o teste estatístico adequado para as comparações (2 grupos)
| Distribuição | Teste | Aplicação |
|---|---|---|
| Paramétrico | Teste t de Student para amostras independentes | Compara as médias de dois grupos independentes quando a distribuição dos dados é normal e as variâncias são iguais. |
| Paramétrico | Teste t de Student para amostras pareadas | Compara as médias de dois grupos pareados ou dependentes quando a distribuição dos dados é normal. |
| Paramétrico | Teste t de Welch | Compara as médias de dois grupos independentes quando a distribuição dos dados é normal, mas as variâncias são diferentes. |
| Não Paramétrico | Teste de Mann-Whitney | Compara as medianas de dois grupos independentes quando a distribuição dos dados não é normal. |
| Não Paramétrico | Teste de Wilcoxon para amostras pareadas | Compara as medianas de dois grupos pareados ou dependentes quando a distribuição dos dados não é normal. |
| Não Paramétrico | Teste de Kolmogorov-Smirnov | Compara as distribuições de dois grupos independentes, sem assumir uma distribuição específica. |
| Não Paramétrico | Teste de Fisher | Compara proporções entre dois grupos independentes quando as amostras são pequenas. |
Desse modo, podemos atribuir a seguinte sequência de testes para as comparações.
resumo |> dplyr::mutate(
"Teste Estatístico" = dplyr::case_when(
`Amostra Feminina` == "Paramétrico" & `Amostra Masculina` == "Paramétrico" ~ "Teste da Família t",
TRUE ~ "Teste de Mann-Whitney ou Wilcoxon")) |>
gt::gt() |>
gt::tab_header(gt::md("**Escolhendo o teste estatístico com base na distribuição dos dados**"),
subtitle = "Prof. Dr. Gleidson Rebouças") |>
gt::opt_align_table_header(align = "left") |>
gt::cols_label(`Amostra Feminina` = "Feminina",
`Amostra Masculina` = "Masculina") |>
gt::tab_spanner(label = "Amostra", columns = c(2,3)) |>
gt::tab_style(style = gt::cell_text(weight = "bold"),
locations = list(gt::cells_column_labels(),
gt::cells_column_spanners())) |>
gt::cols_align(align = "center", columns = 1) |>
gt::tab_options(table.width = "100%")| Escolhendo o teste estatístico com base na distribuição dos dados | |||
| Prof. Dr. Gleidson Rebouças | |||
| Variáveis Dependentes | Amostra | Teste Estatístico | |
|---|---|---|---|
| Feminina | Masculina | ||
| age | Não Paramétrico | Paramétrico | Teste de Mann-Whitney ou Wilcoxon |
| mRR | Paramétrico | Paramétrico | Teste da Família t |
| SDNN | Paramétrico | Não Paramétrico | Teste de Mann-Whitney ou Wilcoxon |
| rMSSD | Não Paramétrico | Não Paramétrico | Teste de Mann-Whitney ou Wilcoxon |
| LF | Não Paramétrico | Não Paramétrico | Teste de Mann-Whitney ou Wilcoxon |
| LFnu | Paramétrico | Paramétrico | Teste da Família t |
| HF | Não Paramétrico | Não Paramétrico | Teste de Mann-Whitney ou Wilcoxon |
| HFnu | Paramétrico | Paramétrico | Teste da Família t |
| LF_HF | Não Paramétrico | Não Paramétrico | Teste de Mann-Whitney ou Wilcoxon |
Com base na tabela acima podemos estabelecer uma série de decisões. Por exemplo, quais os testes que vamos utilizar e se precisaremos de mais alguma verificação que seja exigência ao teste como por exemplo o teste de igualdade de variâncias de Levene.
Já sabemos que vamos utilizar o teste t para amostras independentes em 3 comparações em função do sexo (mRR, LFnu e HFnu). Em todas as demais iremos utilizar o teste de Mann-Whitney que é o equivalente paramétrico para o teste t independente. Então vamos conhecer algumas informações sobre esse teste antes de continuar.
O teste de Mann-Whitney e o teste de Wilcoxon rank-sum são testes não-paramétricos equivalentes que comparam duas amostras independentes para determinar se elas vêm da mesma distribuição. Este último ainda é diferente do teste de Wilconxon (tradicional) que é utilizado para amostras pareadas como alternativa não paramétrica ao teste t pareado. (Lehmann 2011; Hollander, A. Wolfe, e Chicken 2015; Wilcoxon 1945; Mann e Whitney 1947)
Ambos os testes são matematicamente equivalentes e frequentemente usados de forma intercambiável.Teste de Mann-Whitney enfatiza a estatística U, enquanto o teste de Wilcoxon enfatiza a soma dos ranks.
Nota
Teste de Mann-Whitney
- Desenvolvedores: Henry B. Mann e Donald R. Whitney, 1947.
- Outro nome: Teste U de Mann-Whitney.
- Objetivo: Comparar duas amostras independentes.
- Enfoque: Estatística U, baseada nas somas das posições (ranks) das observações nas duas amostras.
Nota
Teste de Wilcoxon Rank-Sum
- Desenvolvedor: Frank Wilcoxon, 1945.
- Objetivo: Comparar duas amostras independentes.
- Enfoque: Soma das posições (ranks) das observações nas amostras.
O pacote gtsummary tem uma função que realiza testes de comparação de acordo com a quantidade de grupos a serem comparados. Contudo, os testes realizados como default nem sempre satisfazem nossas intenções. Vejamos abaixo.
gtsummary::reset_gtsummary_theme()
gtsummary::theme_gtsummary_language("pt")
HRV |> dplyr::select(-phys_Act) |>
gtsummary::tbl_summary(by = sex) |>
gtsummary::add_p(pvalue_fun = ~ gtsummary::style_pvalue(.x, digits = 3)) |>
gtsummary::as_gt() |>
gt::tab_options(table.width = "100%")| Características | F, N = 2391 | M, N = 2611 | Valor p2 |
|---|---|---|---|
| Idade (anos) | 31 (28, 34) | 30 (27, 33) | 0.082 |
| mRR | 936 (863, 984) | 930 (867, 986) | 0.898 |
| SDNN | 50 (39, 60) | 50 (37, 59) | 0.621 |
| rMSSD | 40 (23, 56) | 40 (27, 56) | 0.732 |
| LF | 519 (378, 711) | 506 (336, 712) | 0.526 |
| LFnu | 51 (45, 58) | 52 (46, 59) | 0.262 |
| HF | 554 (116, 1,176) | 640 (83, 1,114) | 0.941 |
| HFnu | 40 (33, 46) | 39 (33, 46) | 0.786 |
| LF/HF | 1.26 (1.26, 2.64) | 1.26 (1.26, 2.34) | 0.465 |
| 1 Mediana (AIQ) | |||
| 2 Teste de soma de postos de Wilcoxon | |||
Repare que todas as comparações foram realizadas com o teste de soma dos postos de Wilcoxon e, confome já falamos antes, existem 3 comparações que gostaríamos de realizar com o teste t de Student (mRR, LFnu e HFnu). Desse modo, temos que implementar um parâmetro para informar ao gtsummary que desejamos alterar o teste default para alguma das comparações. Veja abaixo!
gtsummary::reset_gtsummary_theme()
gtsummary::theme_gtsummary_language("pt")
HRV |> dplyr::select(-phys_Act) |>
gtsummary::tbl_summary(by = sex) |>
gtsummary::add_p(pvalue_fun = ~ gtsummary::style_pvalue(.x, digits = 3),
test = list(mRR ~ "t.test", LFnu ~ "t.test", HFnu ~ "t.test")) |>
gtsummary::as_gt() |>
gt::tab_options(table.width = "100%")| Características | F, N = 2391 | M, N = 2611 | Valor p2 |
|---|---|---|---|
| Idade (anos) | 31 (28, 34) | 30 (27, 33) | 0.082 |
| mRR | 936 (863, 984) | 930 (867, 986) | 0.979 |
| SDNN | 50 (39, 60) | 50 (37, 59) | 0.621 |
| rMSSD | 40 (23, 56) | 40 (27, 56) | 0.732 |
| LF | 519 (378, 711) | 506 (336, 712) | 0.526 |
| LFnu | 51 (45, 58) | 52 (46, 59) | 0.303 |
| HF | 554 (116, 1,176) | 640 (83, 1,114) | 0.941 |
| HFnu | 40 (33, 46) | 39 (33, 46) | 0.799 |
| LF/HF | 1.26 (1.26, 2.64) | 1.26 (1.26, 2.34) | 0.465 |
| 1 Mediana (AIQ) | |||
| 2 Teste de soma de postos de Wilcoxon; Teste t com correção de Welch | |||
Conforme pode ser observado no rodapé, temos agora a execução de 2 testes. Se você comparar as 2 tabelas, verá que somente as 3 variáveis que foram informadas quanto à utilização de um outro teste (teste t) houve mudanças nos valores de p. Perceba ainda que o teste t empregado foi com a correção de Welch. Este teste levou em consideração que não houve igualdade de variâncias nos grupos. O teste de homogeneidade de variância de Levene com o pacote car pode nos revelar se a aplicação da correção foi adequada.
Nota
Hipóteses do Teste de Levene
- H0: As variâncias dos grupos são homogêneas (p > 0,05)
- H1: As variâncias dos grupos não são homogêneas (p < 0,05)
library(car)
car::leveneTest(mRR ~ sex, HRV, center = mean)Levene's Test for Homogeneity of Variance (center = mean)
Df F value Pr(>F)
group 1 1.5886 0.2081
498 car::leveneTest(LFnu ~ sex, HRV, center = mean)Levene's Test for Homogeneity of Variance (center = mean)
Df F value Pr(>F)
group 1 0.2726 0.6018
498 car::leveneTest(HFnu ~ sex, HRV, center = mean)Levene's Test for Homogeneity of Variance (center = mean)
Df F value Pr(>F)
group 1 0.0425 0.8368
498
O teste de Levene revelou que as variâncias foram iguais e, portanto, não era necessário a correção de Welch. Sendo assim, podemos mudar parâmetros dentro da função add_p() que informará que as variâncias são iguais. Observe os parâmetros test = e test.args =.
gtsummary::reset_gtsummary_theme()
gtsummary::theme_gtsummary_language("pt")
gtsummary::theme_gtsummary_mean_sd()
HRV |> dplyr::select(-phys_Act) |>
gtsummary::tbl_summary(by = sex,
digits = list(everything() ~ 2)) |>
gtsummary::add_p(
test = list(
age ~ "wilcox.test",
mRR ~ "t.test",
SDNN ~ "wilcox.test",
rMSSD ~ "wilcox.test",
LF ~ "wilcox.test",
LFnu ~ "t.test",
HF ~ "wilcox.test",
HFnu ~ "t.test",
LF_HF ~ "wilcox.test"),
test.args = list(
mRR = list(var.equal = TRUE),
LFnu = list(var.equal = TRUE),
HFnu = list(var.equal = TRUE)),
pvalue_fun = ~ gtsummary::style_pvalue(.x, digits = 3)) |>
gtsummary::as_gt() |>
gt::tab_options(table.width = "100%")| Características | F, N = 2391 | M, N = 2611 | Valor p2 |
|---|---|---|---|
| Idade (anos) | 30.40 (4.47) | 29.85 (4.60) | 0.082 |
| mRR | 928.63 (93.35) | 928.84 (84.06) | 0.979 |
| SDNN | 51.07 (14.45) | 50.24 (14.03) | 0.621 |
| rMSSD | 42.63 (20.03) | 42.97 (19.40) | 0.732 |
| LF | 546.19 (241.24) | 534.34 (250.97) | 0.526 |
| LFnu | 51.39 (9.74) | 52.29 (9.72) | 0.303 |
| HF | 732.92 (624.89) | 738.13 (603.41) | 0.941 |
| HFnu | 39.67 (10.06) | 39.44 (9.93) | 0.798 |
| LF/HF | 2.73 (2.76) | 2.66 (2.73) | 0.465 |
| 1 Média (Desvio Padrão) | |||
| 2 Teste de soma de postos de Wilcoxon; Teste t para amostras independentes | |||
Note que nas notas de rodapé, não consta mais a correção de Welch, dado o fato de que informamos não haver diferenças nas variâncias.
Observando agora os valores de p para cada uma das comparações, podemos perceber que não encontramos diferenças significativas nas comparações de média para todas as variáveis.
Se por acaso não ficou claro a utilização do teste t dentro dos parâmetros do pacote gtsummary, você pode realizar os testes separadamente usando a função t.test nativa do R e em seguida copiar e alimentar uma tabela à medida que for obtendo as informações. Veja abaixo como fazer:
# Função nativa do R para o teste t de student
t.test(mRR ~ sex, HRV, var.equal = TRUE)
Two Sample t-test
data: mRR by sex
t = -0.02588, df = 498, p-value = 0.9794
alternative hypothesis: true difference in means between group F and group M is not equal to 0
95 percent confidence interval:
-15.79459 15.38389
sample estimates:
mean in group F mean in group M
928.6315 928.8369 t.test(LFnu ~ sex, HRV, var.equal = TRUE)
Two Sample t-test
data: LFnu by sex
t = -1.0319, df = 498, p-value = 0.3026
alternative hypothesis: true difference in means between group F and group M is not equal to 0
95 percent confidence interval:
-2.6095687 0.8123241
sample estimates:
mean in group F mean in group M
51.39180 52.29042 t.test(HFnu ~ sex, HRV, var.equal = TRUE)
Two Sample t-test
data: HFnu by sex
t = 0.25553, df = 498, p-value = 0.7984
alternative hypothesis: true difference in means between group F and group M is not equal to 0
95 percent confidence interval:
-1.528773 1.985883
sample estimates:
mean in group F mean in group M
39.66695 39.43839 # Função nativa do R para o teste t de Wilcoxon
wilcox.test(age ~ sex, HRV)
Wilcoxon rank sum test with continuity correction
data: age by sex
W = 33989, p-value = 0.08216
alternative hypothesis: true location shift is not equal to 0wilcox.test(SDNN ~ sex, HRV)
Wilcoxon rank sum test with continuity correction
data: SDNN by sex
W = 31987, p-value = 0.6212
alternative hypothesis: true location shift is not equal to 0wilcox.test(rMSSD ~ sex, HRV)
Wilcoxon rank sum test with continuity correction
data: rMSSD by sex
W = 30638, p-value = 0.7321
alternative hypothesis: true location shift is not equal to 0wilcox.test(LF ~ sex, HRV)
Wilcoxon rank sum test with continuity correction
data: LF by sex
W = 32212, p-value = 0.5263
alternative hypothesis: true location shift is not equal to 0wilcox.test(HF ~ sex, HRV)
Wilcoxon rank sum test with continuity correction
data: HF by sex
W = 31070, p-value = 0.9406
alternative hypothesis: true location shift is not equal to 0wilcox.test(LF_HF ~ sex, HRV)
Wilcoxon rank sum test with continuity correction
data: LF_HF by sex
W = 32319, p-value = 0.4647
alternative hypothesis: true location shift is not equal to 0Apresentação gráfica dos resultados
Apesar de não termos diferenças estatísticas em nenhuma das comparações, iresmo agora exibir gráficos para melhor visualização disso considerando que, como um tutorial, devemos estar prontos para utilizar de gráficos para uma melhor exibição dos resultados. Para mais informações acesse: https://www.data-to-viz.com/
De acordo com o algorítmo acima, temos algumas possibilidades de exição para as nossas variáveis do dataset HRV. Como já fizemos exibição de linhas de código que envolvem o boxplot anteriormente, iremos fazer a exibição utilizando um gráfico violino (violin plot).
# Criando uma função personalizada
criando_graficos2 <- function(var_name) {
medias_HRV <- HRV |>
dplyr::group_by(sex) |>
dplyr::summarize(media = mean(HRV[[var_name]], na.rm = TRUE))
# !!sym foi necessário pois o argumento precisará ser dinâmico
ggplot(HRV, aes(x = sex, y = !!sym(var_name), fill = sex)) +
geom_violin(trim = FALSE) +
geom_boxplot(width = 0.1) +
scale_fill_brewer(palette = "Greens") +
theme(legend.position = "none")
}
# Selecionar as variáveis numéricas
numeric_vars2 <- names(HRV)[sapply(HRV, is.numeric)]
# Criando os gráficos para todas as variáveis numéricas
meus_plots2 <- purrr::map(numeric_vars2, criando_graficos2)
# Juntando os gráficos
gridExtra::grid.arrange(grobs = meus_plots2, nrow = 3, ncol = 3)
Referências
Bulmer, Michael George. 1979. Principles of statistics. Courier Corporation.
Gupta, Anshul, Prabhaker Mishra, ChandraM Pandey, Uttam Singh, Chinmoy Sahu, e Amit Keshri. 2019. «Descriptive Statistics and Normality Tests for Statistical Data». Annals of Cardiac Anaesthesia 22 (1): 67. https://doi.org/10.4103/aca.aca_157_18.
Hair, Joseph F. 2009. «Multivariate data analysis».
Hollander, Myles, Douglas A. Wolfe, e Eric Chicken. 2015. «Nonparametric Statistical Methods». Wiley Series in Probability and Statistics, julho. https://doi.org/10.1002/9781119196037.
Kline, Rex B. 2023. Principles and practice of structural equation modeling. Guilford publications.
Lehmann, Erich L. 2011. «Parametric Versus Nonparametrics: Two Alternative Methodologies». Em, 437–45. Springer US. https://doi.org/10.1007/978-1-4614-1412-4_38.
Lilliefors, Hubert W. 1967. «On the Kolmogorov-Smirnov Test for Normality with Mean and Variance Unknown». Journal of the American Statistical Association 62 (318): 399–402. https://doi.org/10.1080/01621459.1967.10482916.
Mallery, Paul, e Darren George. 2000. SPSS for windows step by step. Allyn & Bacon, Inc.
Mann, H. B., e D. R. Whitney. 1947. «On a Test of Whether One of Two Random Variables Is Stochastically Larger Than the Other». The Annals of Mathematical Statistics 18 (1): 50–60. https://doi.org/10.1214/aoms/1177730491.
Moore, David S, George P McCabe, e Bruce A Craig. 2009. Introduction to the Practice of Statistics. Vol. 4. WH Freeman New York.
Ptaszyński, Krzysztof. 2022. «Bounds on Skewness and Kurtosis of Steady-State Currents». Physical Review E 106 (2). https://doi.org/10.1103/physreve.106.024119.
Shapiro, S. S., e M. B. Wilk. 1965. «An Analysis of Variance Test for Normality (Complete Samples)». Biometrika 52 (3/4): 591. https://doi.org/10.2307/2333709.
Wilcoxon, Frank. 1945. «Individual Comparisons by Ranking Methods». Biometrics Bulletin 1 (6): 80. https://doi.org/10.2307/3001968.