Checking syntax for Welch's t-test, Cohen's d calculation, and looking for ggplot2 tips

Hi everyone,

I am working on an academic assignment using R and I want to double-check my logic and syntax to ensure I am following best practices.

Part 1: Welch's t-test and Effect Size I am testing the effect of medication status on seizure frequency. Since my Levene's test was significant (assumption of homogeneity of variances violated), I opted for a Welch's t-test and calculated Cohen's d without a pooled standard deviation. Could you please confirm if using var.equal = FALSE in t.test() and pooled_sd = FALSE in effectsize::cohens_d() is the correct and most robust approach for this scenario?

Part 2: ggplot2 I also created a scatter plot with a linear regression line to check the correlation between Height and Weight. The code works, but I would love to hear if there are any tips or better geoms/themes to make this graph look more professional and polished for an academic paper.

Here is a reproducible example with dummy data:

Reading the data file

haifa <- read.csv(file.choose())

Loading required packages

library(car)
library(ggpubr)
library(dplyr)
library(effectsize)
library(ggplot2)
library(afex)
library(performance)
library(emmeans)

install.packages("reprex") # Run this only once if needed

---------------------------------------------------------

Question 1 - Independent Samples t-test

The effect of medication status on seizure frequency

---------------------------------------------------------

Checking assumptions for the variables

1 - Checking the normality assumption

First, let's look at the QQ plots for both groups

ggqqplot(haifa$Seizure.Frequency[haifa$Medication.Status == "On Medication"])
ggqqplot(haifa$Seizure.Frequency[haifa$Medication.Status == "Not on Medication"])

Shapiro-Wilk normality test

shapiro.test(haifa$Seizure.Frequency[haifa$Medication.Status == "On Medication"])
shapiro.test(haifa$Seizure.Frequency[haifa$Medication.Status == "Not on Medication"])

2 - Checking the homogeneity of variances assumption

leveneTest(Seizure.Frequency ~ Medication.Status, data = haifa, center = "median")

Independent samples t-test (Welch's t-test due to unequal variances)

t.test(Seizure.Frequency ~ Medication.Status, data = haifa, var.equal = FALSE)

Descriptive statistics

by(haifa$Seizure.Frequency, haifa$Medication.Status, summary)

---------------------------------------------------------

Research Question 2

Is there a positive correlation between patients' height and weight?

---------------------------------------------------------

Checking assumptions

Drawing a QQ plot for each variable

ggqqplot(haifa$Height)
ggqqplot(haifa$Weight)

Shapiro-Wilk normality test

shapiro.test(haifa$Height)
shapiro.test(haifa$Weight)

Calculating the Pearson correlation

cor.test(haifa$Height, haifa$Weight, alternative = "greater")

---------------------------------------------------------

Question 5 - Creating a scatter plot for the correlation

---------------------------------------------------------

The Plot

ggplot(haifa, aes(x = Height, y = Weight)) +
geom_text(label = " :up_down_arrow: ", size = 2) +
geom_smooth(method = "lm", color = "darkgoldenrod", se = TRUE, fill = "grey70") +
stat_cor(method = "pearson", color = "darkgoldenrod", size = 5) +
labs(title = "Is there a positive correlation between Height and Weight?",
x = "Height",
y = "Weight")

---------------------------------------------------------

Question 6 - Effect Size

---------------------------------------------------------

Calculating Cohen's d effect size (un-pooled SD due to Welch's test)

cohens_d(Seizure.Frequency ~ Medication.Status, data = haifa, pooled_sd = FALSE)
link of the data :
https://www.kaggle.com/datasets/amanik000/epilepsy-disorder-dataset

A comment on the ggplot2 plot: I would drop a lot of the details you have added. I'm a minimalist on graphs.

First of all, I am using data.table notation simply because it is easier for me and I have cleaned up the variable names by running the data set through the clean_names() function is {janitor}

I'd drop the regression line, given that data cloud I doubt it adds anything. Rather than using geom_text(label = " :up_down_arrow: ", size = 2) I'd suggest using geom_point(). You have to much data to make the text in geom_text() be of any informational value. All it does is add clutter, or as Tufte calls it, "chartjunk".

So I'd suggest something like this

ibrary(data.table)
library(ggplot2)
library(ggpubr)
library(patchwork)
library(janitor)

# Load data and convert to data.table -------------------------------------
haifa <- read.csv("/home/jrkane/RJunk/Proj_cjhjjdss/raw_data/Epilepsy_dataset.csv") |> clean_names()
DT <- as.data.table(haifa)


# Plot  ---------------------------------------------

ggplot(haifa, aes(x = height, y = weight)) +   geom_point() +
  stat_cor(method = "pearson", color = "darkgoldenrod", size = 5) +
  labs(title = " Height versus  Weight",
       x = "Height",
       y = "Weight")

or to add a bit of colour and interest perhaps

Plot with Sex ----------------------------------------------------------

ggplot(haifa, aes(x = height, y = weight, colour = gender)) + geom_point() +
stat_cor(method = "pearson", color = "darkgoldenrod", size = 5) +
labs(title = " Height versus Weight by Sex",
x = "Height",
y = "Weight")

1 Like

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.

If you have a query related to it or one of the replies, start a new topic and refer back with a link.