---
title: "2025_Daniels_ortho"
author: "Sergio Uribe"
date: 2025-04-10
date-modified: last-modified
language:
title-block-published: "CREATED"
title-block-modified: "UPDATED"
format:
html:
toc: true
number-sections: true
code-fold: true
code-tools: true
toc-expand: 3
theme: cosmo # optional, or pick another like journal, lumen, etc.
toc-expand: 3
code-fold: true
code-tools: true
editor: visual
execute:
echo: false
cache: true
warning: false
message: false
---
# Packages
```{r}
# Load required libraries with pacman; installs them if not already installed
pacman:: p_load (MASS, # for the robust regression
tidyverse, # tools for data science
visdat, #NAs
hrbrthemes, # nice ggplot themes
janitor, # for data cleaning and tables
here, # for reproducible research
gtsummary, # for tables
explore, # for EDA, check https://rolkra.github.io/explore/
# countrycode, # to normalize country data
tidyplots, # check https://tidyplots.org/
ggeasy, # to use variable labels in ggplot, then easy_labs()
easystats, # check https://easystats.github.io/easystats/
scales,
viridis,
flextable,
ggeffects, # for the polot of fits 6 or 12
officer,
writexl, # to save tables in xls format
tidymodels,
sjPlot,
haven, # for spss data
quantreg, # for the quantile regression
lubridate
)
```
```{r}
theme_set (theme_minimal ())
```
# DATA PREPARATION
```{r}
here:: here ()
```
```{r}
df <- read_spss ("ORTOnew.sav" ) |>
select (- 1 ) # Remove names
```
```{r}
# glimpse(df)
```
## Name changes:
```{r}
df <- df |>
rename (
Gender = Dzimums,
Age = Vecums,
Occupation_type = Nodarb_veids,
Education_level = Izglitības_līmenis,
Household_size = Mājsaimniecībā_dzīvoju,
Doctor_communication = Komunikācijā_ar_ārstu,
Family_informed = Mana_ģimene_informēta,
Family_pressure = Spiediens_no_ģimenes,
Doctor_pressure = Spiediens_no_ārsta,
Financial_tracking = Sekoju_finansem,
Extra_budget_money = Naudu_arpus_budz,
Psych_emotional_symptoms = Psihoemoc_trauc,
BFI_Conscientiousness = BFI_Apzinīgums,
BFI_Fake = BFI_Viltus,
BFI_Neuroticism = BFI_Neirotisms,
Aligns_6mo = Kape_atbislt_6men,
Aligns_12mo = Kape_atbilst_12men,
Fits_6mo = Kape_pieguļ_6men,
Fits_12mo = Kape_pieguļ_12men,
Patient_missed_6mo = Pacients_kavē_6men,
Patient_missed_12mo = Pacients_kavē_12men,
Parents_involved_6mo = Vecāki_iesaistīti_6men,
Parents_involved_12mo = Vecāki_iesaistīti_12men,
Shows_interest_6mo = Izrāda_interesi_6men,
Shows_interest_12mo = Izrāda_interesi_12men,
Poor_hygiene_6mo = Slikta_higiēna_6men,
Poor_hygiene_12mo = Slikta_higiēna_12men,
Uses_rubberbands_6mo = Izmanto_gumijas_6men,
Uses_rubberbands_12mo = Izmanto_gumijas_12men,
Shows_dissatisfaction_6mo = Izrāda_neapm_6men,
Shows_dissatisfaction_12mo = Izrāda_neapm_12men,
IOTN_AH = IOTN_AH, # Already in English
IOTN_DH = IOTN_DH, # Already in English
Starting_position_mm = Starting_position_mm, # Already in English
Expected_position_mm = Expected_position_mm, # Already in English
Real_position_mm = Real_position_mm # Already in English
)
```
```{r}
# visdat::vis_dat(df)
```
## Data cleaning
Convert labelled values:\
```{r}
df <- df %>%
mutate (across (where (is.labelled), as_factor))
```
### Transform variables: Descriptive
```{r}
df <- df %>%
mutate (Gender = fct_drop (na_if (Gender, "Nevēlos atbildēt" )))
```
```{r}
df <- df %>%
mutate (Occupation_type = fct_recode (Occupation_type,
"Student" = "Mācos" ,
"Employed" = "Strādāju" ,
"Unemployed" = "Nestrādāju" ,
"Other" = "Cits"
)) %>%
mutate (Occupation_type = na_if (Occupation_type, "2.1" )) %>%
mutate (Occupation_type = fct_drop (Occupation_type))
```
```{r}
df <- df %>%
mutate (Occupation_type = fct_collapse (Occupation_type,
"Other/Unknown" = c ("Other" , "Unknown" )
)) %>%
mutate (Occupation_type = fct_relevel (Occupation_type,
"Student" , "Employed" , "Unemployed" , "Other/Unknown"
)) %>%
mutate (Occupation_type = fct_drop (Occupation_type))
```
(Why still unknown is shown???)
```{r}
df <- df %>%
mutate (
Occupation_type = as_factor (Occupation_type),
Occupation_type = fct_drop (Occupation_type)
)
```
ahh, is a missing value!
```{r}
df <- df %>%
mutate (
Occupation_type = fct_na_value_to_level (Occupation_type, level = "Other/Unknown" ),
Occupation_type = fct_collapse (Occupation_type,
"Other/Unknown" = c ("Other/Unknown" , "Unknown" )
),
Occupation_type = fct_relevel (Occupation_type,
"Student" , "Employed" , "Unemployed" , "Other/Unknown"
)
)
```
```{r}
df <- df %>%
mutate (Education_level = fct_recode (Education_level,
"Primary education" = "Pamatizglitība" ,
"Secondary education" = "Vidēja izglitība" ,
"Vocational secondary" = "Vidēja speciāla izglitība" ,
"Higher education" = "Augstāka izglitība"
)) %>%
mutate (Education_level = fct_relevel (Education_level,
"Primary education" ,
"Secondary education" ,
"Vocational secondary" ,
"Higher education"
))
```
```{r}
df <- df %>%
mutate (Household_size = fct_collapse (Household_size,
"With parents" = "Ar vecākiem" ,
"With romantic/marital partner" = "Ar romantisko attiecību vai laulību partneri" ,
"Alone" = "viens(-a) pats(-i)" ,
"Other" = c ("Ar draugiem vai paziņam" , "5" )
)) %>%
mutate (Household_size = fct_relevel (Household_size,
"With parents" , "With romantic/marital partner" , "Alone" , "Other"
))
```
```{r}
df <- df |>
mutate (
Household_size = fct_recode (Household_size,
"With partner" = "With romantic/marital partner"
)
)
```
```{r}
df <- df |>
mutate (
Doctor_communication = as.character (Doctor_communication),
Doctor_communication = case_when (
Doctor_communication %in% c ("Daļēji piekrītu" ) ~ "Somewhat agree" ,
Doctor_communication %in% c ("Noteikti piekrītu" ) ~ "Strongly agree" ,
Doctor_communication %in% c ("Neitrāli" , "Noteikti nepiekrītu" , "Daļēji nepiekrītu" ) ~ "Other opinion" ,
TRUE ~ Doctor_communication
),
Doctor_communication = factor (Doctor_communication,
levels = c ("Somewhat agree" , "Strongly agree" , "Other opinion" )
)
)
```
```{r}
df <- df |>
mutate (
Family_informed = as.character (Family_informed),
Family_informed = case_when (
Family_informed == "Noteikti piekrītu" ~ "Strongly agree" ,
Family_informed %in% c ("Daļēji nepiekrītu" , "Neitrāli" , "Daļēji piekrītu" , "Noteikti nepiekrītu" ) ~ "Other opinion" ,
TRUE ~ Family_informed
),
Family_informed = factor (Family_informed,
levels = c ("Strongly agree" , "Other opinion" )
)
)
```
```{r}
df <- df |>
mutate (
Family_pressure = as.character (Family_pressure),
Family_pressure = case_when (
Family_pressure == "Noteikti nepiekrītu" ~ "Strongly disagree" ,
Family_pressure == "Daļēji nepiekrītu" ~ "Somewhat disagree" ,
Family_pressure == "Neitrāli" ~ "Neutral" ,
Family_pressure %in% c ("Daļēji piekrītu" , "Noteikti piekrītu" ) ~ "Other opinion" ,
TRUE ~ Family_pressure
),
Family_pressure = factor (Family_pressure,
levels = c ("Strongly disagree" , "Somewhat disagree" , "Neutral" , "Other opinion" )
)
)
```
```{r}
df <- df |>
mutate (
Doctor_pressure = as.character (Doctor_pressure),
Doctor_pressure = case_when (
Doctor_pressure == "Noteikti nepiekrītu" ~ "Strongly disagree" ,
Doctor_pressure == "Daļēji nepiekrītu" ~ "Somewhat disagree" ,
Doctor_pressure == "Neitrāli" ~ "Neutral" ,
Doctor_pressure %in% c ("Daļēji piekrītu" , "Noteikti piekrītu" ) ~ "Other opinion" ,
TRUE ~ Doctor_pressure
),
Doctor_pressure = factor (Doctor_pressure,
levels = c ("Strongly disagree" , "Somewhat disagree" , "Neutral" , "Other opinion" )
)
)
```
```{r}
df <- df |>
mutate (
Financial_tracking = as.character (Financial_tracking),
Financial_tracking = case_when (
Financial_tracking %in% c ("Noteikti nepiekrītu" , "Daļēji nepiekrītu" ) ~ "Other opinion" ,
Financial_tracking == "Neitrāli" ~ "Neutral" ,
Financial_tracking == "Daļēji piekrītu" ~ "Somewhat agree" ,
Financial_tracking == "Noteikti piekrītu" ~ "Strongly agree" ,
TRUE ~ Financial_tracking
),
Financial_tracking = factor (Financial_tracking,
levels = c ("Strongly agree" , "Somewhat agree" , "Neutral" , "Other opinion" )
)
)
```
```{r}
df <- df |>
mutate (
Extra_budget_money = as.character (Extra_budget_money),
Extra_budget_money = case_when (
Extra_budget_money == "Daļēji piekrītu" ~ "Somewhat agree" ,
TRUE ~ "Other opinion"
),
Extra_budget_money = factor (Extra_budget_money,
levels = c ("Somewhat agree" , "Other opinion" )
)
)
```
```{r}
df <- df |>
mutate (
Psych_emotional_symptoms = case_when (
Psych_emotional_symptoms == "Nē" ~ "No" ,
Psych_emotional_symptoms %in% c ("Jā" , "Nevēlos atbildēt" ) ~ "Yes or skipped" ,
TRUE ~ Psych_emotional_symptoms
),
Psych_emotional_symptoms = factor (Psych_emotional_symptoms,
levels = c ("No" , "Yes or skipped" )
)
)
```
### Transform variables: Doctor perception
```{r}
df <- df |>
mutate (
Aligns_12mo = case_when (
Aligns_12mo == "Jā" ~ "Yes" ,
Aligns_12mo == "Nē" ~ "No" ,
TRUE ~ Aligns_12mo
),
Aligns_12mo = factor (Aligns_12mo, levels = c ("Yes" , "No" ))
)
```
Fits_6mo and Fits_12mo (merge rare levels into “Poor fit”
```{r}
fit_levels_good <- "Kape pieguļ pilnībā visos sekstantos; nav vērojamās atstarpes starp kapes materiālu un zobu virsmām"
fit_levels_ok <- "Kape pieguļ pilnībā/vērojamā neliela atstarpe starp kapi un 1/2 zobiem, kuri atrodas vienā/divos blakussekstantos."
df <- df |>
mutate (
Fits_6mo = case_when (
Fits_6mo == fit_levels_good ~ "Perfect fit" ,
Fits_6mo == fit_levels_ok ~ "Slight gap (1–2 teeth)" ,
TRUE ~ "Poor fit"
),
Fits_12mo = case_when (
Fits_12mo == fit_levels_good ~ "Perfect fit" ,
Fits_12mo == fit_levels_ok ~ "Slight gap (1–2 teeth)" ,
TRUE ~ "Poor fit"
),
Fits_6mo = factor (Fits_6mo, levels = c ("Perfect fit" , "Slight gap (1–2 teeth)" , "Poor fit" )),
Fits_12mo = factor (Fits_12mo, levels = c ("Perfect fit" , "Slight gap (1–2 teeth)" , "Poor fit" ))
)
```
Patient_missed_6mo and \_ 12mo
```{r}
df <- df |>
mutate (
Patient_missed_6mo = case_when (
Patient_missed_6mo == "Nekad" ~ "Never" ,
TRUE ~ "Missed appointments"
),
Patient_missed_12mo = case_when (
Patient_missed_12mo == "Nekad" ~ "Never" ,
TRUE ~ "Missed appointments"
),
Patient_missed_6mo = factor (Patient_missed_6mo, levels = c ("Never" , "Missed appointments" )),
Patient_missed_12mo = factor (Patient_missed_12mo, levels = c ("Never" , "Missed appointments" ))
)
```
Parents_involved_6mo and \_ 12mo
```{r}
df <- df |>
mutate (
Parents_involved_6mo = case_when (
Parents_involved_6mo == "Nekad" ~ "Never" ,
TRUE ~ "Sometimes involved"
),
Parents_involved_12mo = case_when (
Parents_involved_12mo == "Nekad" ~ "Never" ,
TRUE ~ "Sometimes involved"
),
Parents_involved_6mo = factor (Parents_involved_6mo, levels = c ("Never" , "Sometimes involved" )),
Parents_involved_12mo = factor (Parents_involved_12mo, levels = c ("Never" , "Sometimes involved" ))
)
```
Shows_interest_6mo and \_ 12mo
```{r}
df <- df |>
mutate (
Shows_interest_6mo = case_when (
Shows_interest_6mo == "Vienmēr" ~ "Always" ,
Shows_interest_6mo == "Bieži" ~ "Often" ,
TRUE ~ "Sometimes"
),
Shows_interest_12mo = case_when (
Shows_interest_12mo == "Vienmēr" ~ "Always" ,
Shows_interest_12mo == "Bieži" ~ "Often" ,
TRUE ~ "Sometimes"
),
Shows_interest_6mo = factor (Shows_interest_6mo, levels = c ("Always" , "Often" , "Sometimes" )),
Shows_interest_12mo = factor (Shows_interest_12mo, levels = c ("Always" , "Often" , "Sometimes" ))
)
```
Poor_hygiene_6mo and \_ 12mo
```{r}
df <- df |>
mutate (
Poor_hygiene_6mo = case_when (
Poor_hygiene_6mo == "Nekad" ~ "Never" ,
TRUE ~ "Sometimes"
),
Poor_hygiene_12mo = case_when (
Poor_hygiene_12mo == "Nekad" ~ "Never" ,
TRUE ~ "Sometimes"
),
Poor_hygiene_6mo = factor (Poor_hygiene_6mo, levels = c ("Never" , "Sometimes" )),
Poor_hygiene_12mo = factor (Poor_hygiene_12mo, levels = c ("Never" , "Sometimes" ))
)
```
Uses_rubberbands_6mo and \_ 12mo
```{r}
df <- df |>
mutate (
Uses_rubberbands_6mo = case_when (
Uses_rubberbands_6mo == "Pacientam šobrīd nav paredzēta ārstēšana ar elastīgam gumijam" ~ "Not prescribed" ,
Uses_rubberbands_6mo %in% c ("Vienmēr" , "Bieži" ) ~ "Uses regularly" ,
TRUE ~ "Rarely/never"
),
Uses_rubberbands_12mo = case_when (
Uses_rubberbands_12mo == "Pacientam šobrīd nav paredzēta ārstēšana ar elastīgam gumijam" ~ "Not prescribed" ,
Uses_rubberbands_12mo %in% c ("Vienmēr" , "Bieži" ) ~ "Uses regularly" ,
TRUE ~ "Rarely/never"
),
Uses_rubberbands_6mo = factor (Uses_rubberbands_6mo, levels = c ("Not prescribed" , "Uses regularly" , "Rarely/never" )),
Uses_rubberbands_12mo = factor (Uses_rubberbands_12mo, levels = c ("Not prescribed" , "Uses regularly" , "Rarely/never" ))
)
```
Shows_dissatisfaction_6mo and \_ 12mo
```{r}
df <- df |>
mutate (
Shows_dissatisfaction_6mo = case_when (
Shows_dissatisfaction_6mo == "Nekad" ~ "Never" ,
TRUE ~ "Sometimes"
),
Shows_dissatisfaction_12mo = case_when (
Shows_dissatisfaction_12mo == "Nekad" ~ "Never" ,
TRUE ~ "Sometimes"
),
Shows_dissatisfaction_6mo = factor (Shows_dissatisfaction_6mo, levels = c ("Never" , "Sometimes" )),
Shows_dissatisfaction_12mo = factor (Shows_dissatisfaction_12mo, levels = c ("Never" , "Sometimes" ))
)
```
# MODELLING
```{r}
# prepare the data from the 01_data_clean
fit_levels_good <- "Kape pieguļ pilnībā visos sekstantos; nav vērojamās atstarpes starp kapes materiālu un zobu virsmām"
fit_levels_ok <- "Kape pieguļ pilnībā/vērojamā neliela atstarpe starp kapi un 1/2 zobiem, kuri atrodas vienā/divos blakussekstantos."
```
# TABLE 1
```{r}
df <- df |>
rowid_to_column ("id" )
```
```{r}
df |>
dplyr:: select (Gender: BFI_Neuroticism) |>
gtsummary:: tbl_summary ()
```
Save as Table 1
```{r}
read_docx () |>
body_add_flextable (
df |>
select (Gender: BFI_Neuroticism) |>
tbl_summary () |>
as_flex_table ()
) |>
print (target = here ("figures_tables" , "table_1.docx" ))
```
## Motivation
### Table data doctor
```{r}
df |>
select (Aligns_12mo: Shows_dissatisfaction_12mo) |>
gtsummary:: tbl_summary ()
```
## Change in position
```{r}
df |>
mutate (deviation_mm = Real_position_mm - Expected_position_mm) |>
ggplot (aes (x = deviation_mm)) +
geom_histogram (binwidth = 0.2 , fill = "steelblue" , color = "white" ) +
geom_vline (xintercept = 0 , linetype = "dashed" ) +
labs (
title = "Deviation between Real and Expected Position" ,
x = "Deviation (mm)" ,
y = "Number of patients"
)
```
```{r}
df |>
mutate (
actual_movement = Starting_position_mm - Real_position_mm,
expected_movement = Starting_position_mm - Expected_position_mm
) |>
ggplot (aes (x = expected_movement, y = actual_movement,
color = Starting_position_mm)) +
geom_point (alpha = 0.7 ) +
geom_abline (slope = 1 , intercept = 0 , linetype = "dashed" ) +
labs (
title = "Expected vs. Actual Tooth Movement" ,
x = "Expected Movement (mm)" ,
y = "Actual Movement (mm)"
)
```
```{r}
df |>
mutate (
actual_movement = Starting_position_mm - Real_position_mm,
expected_movement = Starting_position_mm - Expected_position_mm
) |>
ggplot (aes (x = expected_movement, y = actual_movement,
color = BFI_Conscientiousness)) +
geom_point (alpha = 0.7 ) +
geom_abline (slope = 1 , intercept = 0 , linetype = "dashed" ) +
labs (
title = "Expected vs. Actual Tooth Movement by BFI_Conscientiousness" ,
x = "Expected Movement (mm)" ,
y = "Actual Movement (mm)" ,
color = "Conscientiousness Score"
) +
theme_minimal ()
```
```{r}
df |>
mutate (
actual_movement = Starting_position_mm - Real_position_mm,
expected_movement = Starting_position_mm - Expected_position_mm
) |>
ggplot (aes (x = expected_movement, y = actual_movement,
color = BFI_Neuroticism)) +
geom_point (alpha = 0.7 ) +
geom_abline (slope = 1 , intercept = 0 , linetype = "dashed" ) +
labs (
title = "Expected vs. Actual Tooth Movement by BFI_Neuroticism" ,
x = "Expected Movement (mm)" ,
y = "Actual Movement (mm)" ,
color = "Conscientiousness Score"
) +
theme_minimal ()
```
# MODELLING
Now the question is: what explain the difference between the expected and the actual movement?
Neuroticism and conscientiousness personality traits were assessed with the validated Big Five Personality Inventory (BFI)
```
BFI_Conscientiousness
```
```
BFI_Fake
```
```
BFI_Neuroticism
```
Patient cooperation was assessed by evaluating clinical fitting of aligners to dental arches and comparison of the planned and achieved upper 1^st^ premolar expansion (%).
```
Fits_6mo
```
```
Fits_12mo
```
```
Patient_missed_6mo:Shows_dissatisfaction_12mo
```
**Index of Orthodontic Treatment Needed**
the higher the worst the orthodontic status of the patient
```
IOTN_AH
```
```
IOTN_DH
```
and the target variable is the difference between the
```
Expected_position_mm
```
and the
```
Real_position_mm
```
Model the difference between expected and actual movement
## Change from the initial to the final position
```{r}
df |>
mutate (id = row_number ()) |>
pivot_longer (
cols = c (Starting_position_mm, Real_position_mm),
names_to = "Stage" ,
values_to = "Position_mm"
) |>
mutate (Stage = factor (Stage, levels = c ("Starting_position_mm" , "Real_position_mm" ))) |>
ggplot (aes (x = Stage, y = Position_mm, group = id, color = Gender)) +
geom_line (alpha = 0.3 , color = "steelblue" ) +
geom_point (size = 1.8 ) +
labs (
title = "Actual Tooth Movement from Start to Real Position" ,
x = "Stage" ,
y = "Position (mm)"
) +
theme_minimal ()
```
## Preparing variables
```{r}
df <- df |>
mutate (
movement_deviation = Real_position_mm - Expected_position_mm
)
```
```{r}
df |>
ggplot (aes (x = movement_deviation)) +
geom_histogram (bins = 10 )
```
Ensure age and gender will be included
```{r}
df <- df |>
mutate (
Sex = as.factor (Gender),
Age = as.numeric (Age)
)
```
## Model 1: Lineal Regression
### For the movement
```{r}
m1 <- lm (
movement_deviation ~
Age +
Sex +
Patient_missed_6mo +
Shows_dissatisfaction_12mo +
IOTN_AH +
IOTN_DH +
BFI_Conscientiousness +
BFI_Neuroticism,
data = df
)
```
```{r}
m1 |>
gtsummary:: tbl_regression ()
```
```{r}
plot_model (m1,
type = "est" ,
show.values = TRUE ,
value.offset = 0.3 ,
title = "Predictors of Tooth Movement Deviation" ,
axis.labels = NULL ,
transform = NULL ) + # no transformation — interpret raw units (mm)
theme_minimal ()
```
```{r}
plot_model (
m1,
type = "est" ,
show.values = TRUE ,
value.offset = 0.3 ,
sort.est = TRUE ,
title = "Predictors of Tooth Movement Deviation" ,
axis.labels = c (
"Neuroticism" ,
"Conscientiousness" ,
"IOTN DH" ,
"IOTN AH" ,
"Dissatisfaction at 12 mo [Sometimes]" ,
"Missed appointments at 6 mo" ,
"Sex [Female]" ,
"Age"
),
transform = NULL # Keep in original scale (mm deviation)
) +
theme_minimal () +
theme (
plot.title = element_text (size = 14 , face = "bold" ),
axis.title.x = element_text (size = 12 ),
axis.text = element_text (size = 10 )
) +
labs (x = "Regression Coefficient (mm deviation)" , y = NULL )
```
### Explore BFI_Conscientiousness & BFI_Neuroticism
```{r}
df |>
ggplot (aes (x = BFI_Neuroticism, y = movement_deviation)) +
geom_point (alpha = 0.3 ) +
geom_smooth (method = "loess" , se = T, color = "steelblue" ) +
# geom_smooth(method = "lm", se = TRUE, fill = "lightblue", color = "steelblue") +
# geom_quantile(quantiles = 0.5, color = "steelblue") +
# geom_smooth(method = "gam", formula = y ~ s(x), se = FALSE, color = "purple") +
# geom_smooth() +
labs (
x = "Neuroticism Score" ,
y = "Difference in Expected Movement (mm)" ,
title = "Tooth Movement Deviation by Neuroticism Level"
) +
theme_minimal ()
```
```{r}
df |>
ggplot (aes (x = BFI_Conscientiousness, y = movement_deviation)) +
geom_point (alpha = 0.3 ) +
geom_smooth (method = "loess" , se = T, color = "steelblue" ) +
labs (
x = "Conscientiousness Score" ,
y = "Difference in Expected Movement (mm)" ,
title = "Tooth Movement Deviation by Conscientiousness Level"
) +
theme_minimal ()
```
```{r}
check_model (m1)
```
There are some issues:\
- Observed vs Model-Predicted Density: check, maybe robust regression, transforming the outcome, or checking outliers?
- check cases 24, 30, 38, 47)
- Consider robust regression or bootstrapping if inference is sensitive.
### Check outliers
```{r}
# df |> slice(c(24, 30, 38, 47))
```
## Multinomial models
```{r}
pacman:: p_load (nnet)
```
### Model for Fits 6 and
```{r}
df$ Fits_6mo <- factor (df$ Fits_6mo,
levels = c ("Perfect fit" , "Slight gap (1–2 teeth)" , "Poor fit" ))
df$ Sex <- as.factor (df$ Sex)
df$ Patient_missed_6mo <- as.factor (df$ Patient_missed_6mo)
df$ Shows_dissatisfaction_12mo <- as.factor (df$ Shows_dissatisfaction_12mo)
```
```{r}
m_fit6 <- multinom (
Fits_6mo ~
Age +
Sex +
Patient_missed_6mo +
Shows_dissatisfaction_12mo +
IOTN_AH +
IOTN_DH +
BFI_Conscientiousness +
BFI_Neuroticism,
data = df
)
```
```{r}
m_fit6 |>
gtsummary:: tbl_regression ()
```
```{r}
plot_model (m_fit6, type = "est" , sort.est = TRUE )
```
```{r}
# m_fit6 |>
# tbl_regression(
# exponentiate = TRUE,
# label = list(
# Age ~ "Age",
# Sex ~ "Sex",
# Patient_missed_6mo ~ "Missed appointments",
# Shows_dissatisfaction_12mo ~ "Dissatisfaction (12mo)",
# IOTN_AH ~ "IOTN-AH",
# IOTN_DH ~ "IOTN-DH",
# BFI_Conscientiousness ~ "Conscientiousness",
# BFI_Neuroticism ~ "Neuroticism"
# )
# ) |>
# modify_header(label ~ "**Characteristic**") |>
# modify_caption("**Multinomial logistic regression predicting aligner fit at 6 months**")
```
```{r}
tbl_regression (m_fit6, exponentiate = TRUE ) |>
as_flex_table () |> # or as_gt()
flextable:: autofit ()
```
```{r}
# RUN ONLY ONCE
tbl_regression (m_fit6, exponentiate = TRUE ) |>
as_tibble () |>
write_xlsx (path = here ("figures_tables" , "fit6_model.xlsx" ))
```
### Model for Fit 12 months
```{r}
df$ Fits_12mo <- factor (df$ Fits_12mo,
levels = c ("Perfect fit" , "Slight gap (1–2 teeth)" , "Poor fit" )) # check order
# Ccheck predictors
df$ Sex <- as.factor (df$ Sex)
df$ Patient_missed_6mo <- as.factor (df$ Patient_missed_6mo)
df$ Shows_dissatisfaction_12mo <- as.factor (df$ Shows_dissatisfaction_12mo)
```
```{r}
m_fit12 <- multinom (
Fits_12mo ~
Age +
Sex +
Patient_missed_6mo +
Shows_dissatisfaction_12mo +
IOTN_AH +
IOTN_DH +
BFI_Conscientiousness +
BFI_Neuroticism,
data = df
)
```
```{r}
m_fit12 |>
gtsummary:: tbl_regression ()
```
```{r}
tbl_regression (m_fit12 , exponentiate = TRUE ) |>
as_tibble () |>
write_xlsx (path = here ("figures_tables" , "fit12_model.xlsx" ))
```
## Model 2 Robust regression
```{r}
m_robust <- rlm (
movement_deviation ~ Age + Sex +
BFI_Conscientiousness + BFI_Neuroticism +
Patient_missed_6mo + Shows_dissatisfaction_12mo +
IOTN_AH + IOTN_DH,
data = df
)
```
```{r}
m_robust |>
gtsummary:: tbl_regression ()
```
```{r}
plot_model (m_robust,
type = "est" ,
show.values = TRUE ,
sort.est = TRUE ,
axis.labels = c ("Age" , "Sex" , "Conscientiousness" , "Neuroticism" ,
"Missed Appointments" , "Dissatisfaction" ,
"IOTN-AH" , "IOTN-DH" ),
title = "Robust Regression Coefficients" )
```
## Tidymodel
```{r}
model_recipe <- recipe (movement_deviation ~ Age + Sex +
BFI_Conscientiousness + BFI_Neuroticism +
# Fits_6mo + Fits_12mo +
Patient_missed_6mo + Shows_dissatisfaction_12mo +
IOTN_AH + IOTN_DH,
data = df) |>
step_dummy (all_nominal_predictors ()) |>
step_normalize (all_numeric_predictors ())
```
```{r}
# the model
lin_mod <- linear_reg () |>
set_engine ("lm" ) |>
set_mode ("regression" )
```
```{r}
# the workflow
movement_wf <- workflow () |>
add_model (lin_mod) |>
add_recipe (model_recipe)
```
```{r}
# fit the model
movement_fit <- movement_wf |>
fit (data = df)
```
```{r}
movement_fit |>
tidy () |>
arrange (p.value) |>
knitr:: kable (
digits = 3 ,
caption = "Tidy Output of Ordinal Regression Model (Sorted by p-value)"
)
```
```{r}
movement_fit |>
tidy () |>
filter (term != "(Intercept)" ) |>
mutate (
term = recode (term,
"BFI_Neuroticism" = "Neuroticism" ,
"BFI_Conscientiousness" = "Conscientiousness" ,
"IOTN_DH" = "IOTN Dental Health" ,
"IOTN_AH" = "IOTN Aesthetic" ,
"Age" = "Age" ,
"Sex_Sieviešu" = "Sex (Female)" ,
"Patient_missed_6mo_Missed.appointments" = "Missed appointments (6 mo)" ,
"Shows_dissatisfaction_12mo_Sometimes" = "Dissatisfaction (12 mo)"
)
) |>
ggplot (aes (x = reorder (term, estimate), y = estimate)) +
geom_point (size = 3 ) +
geom_errorbar (aes (
ymin = estimate - std.error * 1.5 ,
ymax = estimate + std.error * 1.5
), width = 0.2 ) +
coord_flip () +
labs (
title = "Predictors of Tooth Movement Deviation (95% Confidence Interval)" ,
x = "Predictor" ,
y = "Estimate (mm)"
) +
geom_hline (yintercept = 0 , linetype = "dashed" , color = "gray40" ) +
theme_minimal ()
```
```{r}
ggsave (
filename = here ("figures_tables" , "predictors_tooth_movement.png" ),
plot = last_plot (),
width = 10 ,
height = 5 ,
dpi = 300 ,
bg = "transparent"
)
```
Presentation: <https://docs.google.com/presentation/d/1VqYmsV5PSVvMKO3lzGamVfkDsTMCwQqNK47oTP_s-lw/edit#slide=id.p10>
## Effect on subjective measurements
### Fits 6 months
```{r}
df$ Fits_6mo <- factor (
df$ Fits_6mo,
levels = c ("Perfect fit" , "Slight gap (1–2 teeth)" , "Poor fit" ),
ordered = TRUE
)
```
```{r}
# Model 1: BFI_Neuroticism + Gender
model_neurotic_ord <- polr (Fits_6mo ~ BFI_Neuroticism + Gender, data = df, Hess = TRUE )
```
```{r}
model_neurotic_ord |>
gtsummary:: tbl_regression (exponentiate = T)
```
```{r}
# Predicted probabilities from ordinal logistic model
pred_neurotic <- ggeffects:: ggpredict (model_neurotic_ord, terms = "BFI_Neuroticism" )
# Plot
plot (pred_neurotic) +
labs (
title = "Predicted Fit Quality by Neuroticism Score" ,
x = "BFI Neuroticism Score" ,
y = "Predicted Probability"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Neuroticism_group = cut (
BFI_Neuroticism,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Neuroticism_group, fill = Fits_6mo)) +
geom_bar (position = "stack" ) +
labs (
title = "Fit Quality Distribution by Neuroticism Group" ,
x = "Neuroticism Level" ,
y = "Count" ,
fill = "Fit Quality (6 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Neuroticism_group = cut (
BFI_Neuroticism,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Neuroticism_group, fill = Fits_6mo)) +
geom_bar (position = "fill" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
scale_y_continuous (labels = scales:: percent_format ()) +
labs (
title = "Fit Quality (6 months) Distribution by Neuroticism Group" ,
x = "Neuroticism Level" ,
y = "Proportion" ,
fill = "Fit Quality (6 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Neuroticism_group = cut (
BFI_Neuroticism,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Neuroticism_group, fill = Fits_12mo)) +
geom_bar (position = "stack" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
# scale_y_continuous(labels = scales::percent_format()) +
labs (
title = "Fit Quality (12 months) Distribution by Neuroticism Group" ,
x = "Neuroticism Level" ,
y = "Count" ,
fill = "Fit Quality (12 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Neuroticism_group = cut (
BFI_Neuroticism,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Neuroticism_group, fill = Fits_12mo)) +
geom_bar (position = "fill" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
scale_y_continuous (labels = scales:: percent_format ()) +
labs (
title = "Fit Quality (12 months) Distribution by Neuroticism Group" ,
x = "Neuroticism Level" ,
y = "Proportion" ,
fill = "Fit Quality (12 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Conscientiousness_group = cut (
BFI_Conscientiousness,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Conscientiousness_group, fill = Fits_6mo)) +
geom_bar (position = "stack" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
#scale_y_continuous(labels = scales::percent_format()) +
labs (
title = "Fit Quality (6 months) Distribution by Conscientiousness Group" ,
x = "Conscientiousness_group Level" ,
y = "Count" ,
fill = "Fit Quality (6 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Conscientiousness_group = cut (
BFI_Conscientiousness,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Conscientiousness_group, fill = Fits_6mo)) +
geom_bar (position = "fill" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
scale_y_continuous (labels = scales:: percent_format ()) +
labs (
title = "Fit Quality (6 months) Distribution by Conscientiousness Group" ,
x = "Conscientiousness_group Level" ,
y = "Proportion" ,
fill = "Fit Quality (6 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Conscientiousness_group = cut (
BFI_Conscientiousness,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Conscientiousness_group, fill = Fits_12mo)) +
geom_bar (position = "stack" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
# scale_y_continuous(labels = scales::percent_format()) +
labs (
title = "Fit Quality (12 months) Distribution by Conscientiousness Group" ,
x = "Conscientiousness_group Level" ,
y = "Count" ,
fill = "Fit Quality (12 months)"
) +
theme_minimal ()
```
```{r}
df |>
mutate (Conscientiousness_group = cut (
BFI_Conscientiousness,
breaks = c (0 , 40 , 50 , Inf ),
labels = c ("Low" , "Moderate" , "High" ))
) |>
ggplot (aes (x = Conscientiousness_group, fill = Fits_12mo)) +
geom_bar (position = "fill" ) +
scale_fill_viridis_d (option = "D" ) + # You can try options: "A", "B", "C", "D", "E", "F"
scale_y_continuous (labels = scales:: percent_format ()) +
labs (
title = "Fit Quality (12 months) Distribution by Conscientiousness Group" ,
x = "Conscientiousness_group Level" ,
y = "Proportion" ,
fill = "Fit Quality (12 months)"
) +
theme_minimal ()
```