ML_IV_Independent_Analysis_Part_1_F.Rmd

---
title: "ML_IV_Independent_Analysis_Part_1"
output: html_document
date: "2026-02-20"
---

```{r}
knitr::opts_chunk$set(echo = TRUE)
library(rstatix)
library(tidyverse)
library(tidyr)
library(pastecs)
library(knitr)
library(epitools)
library(Epi)
library(readr) 
library(dplyr)
library(gtsummary)
library(reportRmd)
library(naniar)
library(visdat)
library(mice)
options(scipen=999) 
digits <- 2  
```


```{r}
# Loading the dataset

data1 <-read_csv("ML_Immunization_BL_data_20_Feb_26.csv")

```

### Data wrangling

```{r}
# Examining dataset

glimpse(data1)
```


```{r}
# Checking for missing values
missing_table <- miss_var_summary(data1)
missing_table
```

We have 33 variables with missing values ranging from 12 to 88%. although most of these missing values are due to the nature of the data where some questions were not applicable to all respondents I have imputed the missing values for the purpose of this assignment using MICE.
```{r}
# Imputing missing values

mice_all <- mice(data1, method = "cart", m = 1, maxit = 1) 
mice_all$method

mice_immunz_imp <- complete(mice_all)
write_csv(mice_immunz_imp, "mice_immunz_imp.csv")
```

```{r}
data <- read_csv("mice_immunz_imp.csv")
```


```{r}
# Checking for missing values after imputation
missing_table <- miss_var_summary(data)
missing_table
```

Perfect for now

### Numeric Variables

```{r}
# Summary statistics for numeric variables 

rm_covsum(data=data,
covs=c('Age_Child', 'Age_mother', 'Age_father', 'HH_size', 'Child_under_18', 'Child_under_5', 'HH_income', 'Time_to_imm_site', 'No_ANC', 'No_PNC', 'No_sessions_comp_imm', 'Age_comp_imm'))

```

 
As it can see from the above table, the data needs cleaning. For instance, among others, the maximum age of the father is 514 years, the maximum HH size is 2000, the maximum time required to reach to the immunization site is 500 minutes, and the maximum number of sessions to complete immunization is 2000, etc, which are not meaningful values. Therefore, data cleaning is conducted to handle this. 

```{r}
# Data cleaning for numeric variables: outliers and non-meaningful values. 

data <- data %>%
          mutate(Age_father = case_when(
            Age_father > 76 ~ 76,
            Age_father < 19 ~ 19,
            TRUE ~ Age_father
        ))

data <- data %>%
          mutate(HH_size = case_when(
          HH_size > 15 ~ 15,
          TRUE ~ HH_size
          ))

data <- data %>%
          mutate(Child_under_18 = case_when(
          Child_under_18 < 1 ~ 1,
          Child_under_18 > 10 ~ 10,
          TRUE ~ Child_under_18
          ))

data <- data %>%
          mutate(Child_under_5 = case_when(
          Child_under_5 < 1 ~ 1,
          Child_under_5 > 5 ~ 5,
          TRUE ~ Child_under_5
          ))

data <- data %>%
          mutate(HH_income = case_when(
            HH_income > 50000 ~ 50000,
            HH_income < 500 ~ 500,
            TRUE ~ HH_income
        ))

data <- data %>%
          mutate(Time_to_imm_site = case_when(
            Time_to_imm_site > 120 ~ 120,
            Time_to_imm_site < 5 ~ 5,
            TRUE ~ Time_to_imm_site
        ))

data <- data %>%
          mutate(No_sessions_comp_imm = case_when(
            No_sessions_comp_imm > 12 ~ 12,
            No_sessions_comp_imm < 1 ~ 1,
            TRUE ~ No_sessions_comp_imm
        ))

data <- data %>%
          mutate(Age_comp_imm = case_when(
            Age_comp_imm > 60 ~ 60,
            Age_comp_imm < 1 ~ 1,
            TRUE ~ Age_comp_imm
        ))

```


```{r}
## summary statistics after data cleaning 

rm_covsum(data=data,
covs=c('Age_Child', 'Age_mother', 'Age_father', 'HH_size', 'Child_under_18', 'Child_under_5', 'HH_income', 'Time_to_imm_site', 'No_ANC', 'No_PNC', 'No_sessions_comp_imm', 'Age_comp_imm'))

```


### Other variables

```{r}
# Child sex
data |> count(Sex_Child, sort = TRUE)
```

```{r}
# Recoding
data <- data |>
  mutate(Sex_Child = case_when(
    Sex_Child == 1 ~ "Male",
    Sex_Child == 2 ~ "Female"
  ))
data$Sex_Child <- as.factor(data$Sex_Child)
data |> count(Sex_Child, sort = TRUE)
```


```{r}
# Marital_status 
data |> count(Marital_status, sort = TRUE)
```


```{r}
# Recoding 
data <- data %>%
  mutate(Marital_status = case_when(
    Marital_status == 1 ~ "Married",
    Marital_status == 2 ~ "Not married",
    Marital_status == 3 ~ "Living together",
    Marital_status == 4 ~ "Separated",
    Marital_status == 5 ~ "Divorced",
    Marital_status == 6 ~ "Widowed"
  ))
data$Marital_status <- as.factor(data$Marital_status)
data |> count(Marital_status, sort = TRUE)
```


```{r}
# Literacy status of the mother
data |> count(Mother_read_or_write, sort = TRUE)
```


```{r}
# Recoding
data <- data |>
  mutate(Mother_read_or_write = case_when(
    Mother_read_or_write == 1 ~ "Yes",
    Mother_read_or_write == 0 ~ "No"
  ))
data$Mother_read_or_write <- as.factor(data$Mother_read_or_write)
data |> count(Mother_read_or_write, sort = TRUE)
```


```{r}
# Literacy status of the father
data |> count(Father_read_or_write, sort = TRUE)
```


```{r}
# Recoding
data <- data |>
  mutate(Father_read_or_write = case_when(
    Father_read_or_write == 1 ~ "Yes",
    Father_read_or_write == 0 ~ "No"
  ))
data$Father_read_or_write <- as.factor(data$Father_read_or_write)
data |> count(Father_read_or_write, sort = TRUE)
```


```{r}
# Residence
data |> count(Residence, sort = TRUE)
```


```{r}
# Recoding
data <- data |>
  mutate(Residence = case_when(
    Residence == 1 ~ "Urban",
    Residence == 2 ~ "Rural"
  ))
data$Residence <- as.factor(data$Residence)
data |> count(Residence, sort = TRUE)
```

```{r}
# ANC follow up during pregnancy
data |> count(ANC_Px, sort = TRUE)
```


```{r}
# Recoding
data <- data |>
  mutate(ANC_Px = case_when(
    ANC_Px == 1 ~ "Yes",
    ANC_Px == 0 ~ "No"
  ))
data$ANC_Px <- as.factor(data$ANC_Px)
data |> count(ANC_Px, sort = TRUE)
```

```{r}
# Place of child birth
data |> count(Place_birth, sort = TRUE)
```

```{r}
# Recoding
data <- data |>
  mutate(Place_birth = case_when(
    Place_birth == 1 ~ "Health facility",
    Place_birth == 2 ~ "Home"
  ))
data$Place_birth <- as.factor(data$Place_birth)
data |> count(Place_birth, sort = TRUE)
```


```{r}
# PNC follow up after delivery
data |> count(PNC, sort = TRUE)
```


```{r}
# Recoding
data <- data |>
  mutate(PNC = case_when(
    PNC == 1 ~ "Yes",
    PNC == 0 ~ "No"
  ))
data$PNC <- as.factor(data$PNC)
data |> count(PNC, sort = TRUE)
```

```{r}
# Receiving health education about immunization
data |> count(HE_Imm, sort = TRUE)
```

```{r}
# Recoding
data <- data |>
  mutate(HE_Imm = case_when(
    HE_Imm == 1 ~ "Yes",
    HE_Imm == 0 ~ "No"
  ))
data$HE_Imm <- as.factor(data$HE_Imm)
data |> count(HE_Imm, sort = TRUE)
```

```{r}
# When the child should start immunization
data |> count(When_to_immun, sort = TRUE)
```


```{r}
# Recoding
data <- data |>
  mutate(When_to_immun = case_when(
    When_to_immun == 1 ~ "At birth",
    When_to_immun == 2 ~ "After 45 days",
    When_to_immun == 3 ~ "First two weeks",
    When_to_immun == 4 ~ "After two weeks",
    When_to_immun == 5 ~ "Don't know"
  ))
data $When_to_immun <- as.factor(data$When_to_immun)
data |> count(When_to_immun, sort = TRUE)
```

**Note:** "When to immun" is from those who received health education about immunization- Missing is expected here. 

```{r}
# Comprehensive table with all the above variables

rm_covsum(
  data = data,
  covs = c('Sex_Child', 'Marital_status', 'Mother_read_or_write', 'Father_read_or_write', 'Residence', 'ANC_Px', 'Place_birth', 'PNC', 'HE_Imm', 'When_to_immun'),
  show.tests = TRUE
)
```


Most of the below variables are factors with 1 representing Yes and 0 represents No. These variables are changed to factors and presented accordingly. 

```{r}
# Knowledge about vaccine preventable diseases

## Converting to factor

data <- data %>%
  mutate(across(
    c(Identify_Ds_Tb, Identify_Ds_polio, Identify_Ds_diptheria, 
      Identify_Ds_purtusis, Identify_Ds_tetanus, Identify_Ds_measeles, 
      Identify_Ds_pneumonia, Identify_Ds_Meng, Identify_Ds_HPB, 
      Identify_Ds_diarreah),
    ~ factor(.)
  ))
```

```{r}
## Table with all variables relevant to knowledge about vaccine preventable diseases

rm_covsum(
  data = data,
  covs = c('Identify_Ds_Tb', 'Identify_Ds_polio', 'Identify_Ds_diptheria', 'Identify_Ds_purtusis', 'Identify_Ds_tetanus', 'Identify_Ds_measeles', 'Identify_Ds_pneumonia', 'Identify_Ds_Meng', 'Identify_Ds_HPB', 'Identify_Ds_diarreah'), 
  show.tests = TRUE
)
```


```{r}
# Perceptions about immunization: Who should not be immunized_Sick children, Physically handicapped children

## Converting to factor
data <- data |>
  mutate(across(
    c(Some_children_should_not, Who_Newborns, Who_Sick, Who_Physically_handi_ca),
    ~ factor(.)
  ))
```

```{r}
## Table for perceptions about immunization
rm_covsum(
  data = data,
  covs = c('Some_children_should_not', 'Who_Newborns', 'Who_Sick', 'Who_Physically_handi_ca'), 
  show.tests = TRUE
)
```

```{r}
# perceptions about immunization: Vaccines are harmful, Child fever after vaccination, Child diarrhea after vaccination, Doubt about vaccine, Worry about vaccine causing illness, Recommending others to vaccinate

## Converting to factor
data <- data |>
  mutate(across(
    c(Vacc_harmful, Child_fever_vacc, Child_diarrhea_Vacc, Doubt_vacc, Worry_Vacc_illness, Recommending_others),
    ~ factor(.)
  ))
```

```{r}
## Table for perceptions about immunization: 
rm_covsum(
  data = data,
  covs = c('Vacc_harmful', 'Child_fever_vacc', 'Child_diarrhea_Vacc', 'Doubt_vacc', 'Worry_Vacc_illness', 'Recommending_others'), 
  show.tests = TRUE
)
```

```{r}
# Motivation to vaccinate

## Converting to factor
data <- data |>
  mutate(across(
    c(Motivation_prevents_dis, Motivation_prevents_death, Motivation_HW_told),
    ~ factor(.)
    ))
```

```{r}
## Table for motivation to vaccinate
rm_covsum(
  data = data,
  covs = c('Motivation_prevents_dis', 'Motivation_prevents_death', 'Motivation_HW_told'), 
  show.tests = TRUE
)
```


```{r}
# Vaccination status of the child

## Converting to factor
data <- data |>
  mutate(across(
    c(Child_ever_vaccinated, Scar_BCG, vaccination_card, see_card),
    ~ factor(.)
  ))
```

```{r}
## Table for vaccination status of the child
rm_covsum(
  data = data,
  covs = c('Child_ever_vaccinated', 'Scar_BCG', 'vaccination_card', 'see_card'), 
  show.tests = TRUE
)
```


```{r}
# Received specific vaccines 

## Converting to factor
data <- data |>
  mutate(across(
    c(Card_BCG, Card_Polio0, Card_Polio1, Card_Polio2, Card_Polio3, Card_Penta1, Card_Penta2, Card_Penta3, Card_PCV1, Card_PCV2, Card_PCV3, Card_Rota1,  Card_Rota2, Card_MCV1, Card_MCV2),
    ~ factor(.)
  ))
```

```{r}
## Table for receiving specific vaccines
rm_covsum(
  data = data,
  covs = c('Card_BCG', 'Card_Polio0', 'Card_Polio1', 'Card_Polio2', 'Card_Polio3', 'Card_Penta1', 'Card_Penta2', 'Card_Penta3', 'Card_PCV1', 'Card_PCV2', 'Card_PCV3', 'Card_Rota1',  'Card_Rota2', 'Card_MCV1', 'Card_MCV2'),
  show.tests = TRUE
)
```

Problem in MCV2 and MCV 2 where we see a value of 2 which is not meaningful. Therefore, we will recode 2 into 0 as 1 and 0 or 1 and 2 are usually used to represent yes and no in the survey in the area.


```{r}
# Recoding MCV1 and MCV2

data <- data |>
  mutate(Card_MCV1 = case_when(
    Card_MCV1 == "2" ~ "0",
    TRUE ~ as.character(Card_MCV1)
  ),
Card_MCV1 = factor(Card_MCV1, levels = c("0", "1"))
  )

data <- data %>%
  mutate(
    Card_MCV2 = case_when(
      Card_MCV2 == "2" ~ "0",
      TRUE ~ as.character(Card_MCV2)
    ),
    Card_MCV2 = factor(Card_MCV2, levels = c("0", "1"))
  )
```

```{r}
## Check- table for receiving specific vaccines
rm_covsum(
  data = data,
  covs = c('Card_BCG', 'Card_Polio0', 'Card_Polio1', 'Card_Polio2', 'Card_Polio3', 'Card_Penta1', 'Card_Penta2', 'Card_Penta3', 'Card_PCV1', 'Card_PCV2', 'Card_PCV3', 'Card_Rota1',  'Card_Rota2', 'Card_MCV1', 'Card_MCV2'),
  show.tests = TRUE
)
```


#### Dependent variable  

At this stage, I have prepared the following three dependent variables from the dataset. The choice will made in the later stage- Independent data analysis Part 2. 


**1- Continuous Dependent variable:** It is created by adding all the vaccinations the child has received excluding polio 0 and MCV 2 following the recommended practice for the age group of 12-23 months. 

```{r}
## Adding all the vaccinations the child has received excluding polio 0 and MCV 2

vax_cols <- c(
  "Card_BCG", "Card_Polio1", "Card_Polio2", "Card_Polio3",
  "Card_Penta1", "Card_Penta2", "Card_Penta3",
  "Card_PCV1", "Card_PCV2", "Card_PCV3",
  "Card_Rota1", "Card_Rota2", "Card_MCV1"
)

data <- data %>%
  mutate(
    across(all_of(vax_cols), ~ as.numeric(as.character(.))),
    vaccine_count = rowSums(across(all_of(vax_cols)), na.rm = TRUE)
  )
summary(data$vaccine_count)

```
**2- Binary Dependent variable:** Created by categorizing the number of vaccines received into two categories: 

* Fully immunized: Received all the recommended 13 doses of vaccines
* Not fully-immunized: Received at least one dose of the recommended vaccines 

```{r}
# Recoding into two categories
data <- data |>
  mutate(vaccine_count_cat = case_when(
  vaccine_count == 13 ~ "Fully immunized",
  vaccine_count < 13 ~ "Not fully-immunized"
  ))
data$vaccine_count_cat <- as.factor(data$vaccine_count_cat)
```


```{r}
data |>
  count(vaccine_count_cat, .drop = FALSE) |>
  mutate(
    pct = round(100 * n / sum(n), digits)  
  )
```

**3- Binary Dependent variable 2 (Zero doze):** Children who did not receive Penta 1, which is a proxy for not receiving any dose of the recommended vaccines.

```{r}
# Recoding zero doze variable
data <- data |>
  mutate(zero_dose = case_when(
    Card_Penta1 == "0" ~ "Zero doze",
    Card_Penta1 == "1" ~ "Not zero doze"
  ))
data$zero_dose <- as.factor(data$zero_dose)
```

```{r}
data |>
  count(zero_dose, .drop = FALSE) |>
  mutate(
    pct = round(100 * n / sum(n), digits)  
  )
```

Based on the table above, only 5.2% of children are fully immunized. This low coverage is expected because the data was collected in some of the most remote areas of the country and included children from refugee communities. In addition, immunization status in this dataset is based solely on information recorded on vaccination cards, which rarely available in rural households and refugee settings in Ethiopia. The same for Zero doze.




```{r}
# The place to get vaccination
## Converting to factor
data <- data |>
  mutate(across(
    c(Where_outreach, Where_HP, Where_HC, Where_Hospital, Where_Private, Where_Faith_clinic, Where_Other),
    ~ factor(.)
  ))

```

```{r}
## Table for the place to vaccinate the child
rm_covsum(
  data = data,
  covs = c('Where_outreach', 'Where_HP', 'Where_HC', 'Where_Hospital', 'Where_Private', 'Where_Faith_clinic', 'Where_Other'), 
  show.tests = TRUE
)
```



```{r}
# Information given during vaccination service
## Converting to factor
data <- data |>
  mutate(across(
    c(Told_adverse_events, Told_what_to_do, Told_next_visit),
    ~ factor(.)
  ))
```

```{r}
## Table for information given during vaccination service
rm_covsum(
  data = data,
  covs = c('Told_adverse_events', 'Told_what_to_do', 'Told_next_visit'), 
  show.tests = TRUE
)
```


```{r}
# Experience in the vaccination service
## Converting to factor
data <- data |>
  mutate(across(
    c(Returned_without_vacc, Other_HF_visit, Received_vacc, Refused_vaccinations),
    ~ factor(.)
  ))
```

```{r}
## Table for experience in the vaccination service
rm_covsum(
  data = data,
  covs = c('Returned_without_vacc', 'Other_HF_visit', 'Received_vacc', 'Refused_vaccinations'),
  show.tests = TRUE
)
```


```{r}
# Service satisfaction
data |> count(Service_satisfaction, sort = TRUE)
```

```{r}
# Recoding 1 and 2 into satisfied, 3 into neutral and 4 and 5 into dissatisfied
data <- data |>
  mutate(Service_satisfaction = case_when(
    Service_satisfaction %in% c(1, 2) ~ 1,
    Service_satisfaction == 3 ~ 2,
    Service_satisfaction %in% c(4, 5) ~ 3
  ))
data$Service_satisfaction <- as.factor(data$Service_satisfaction)
data |> count(Service_satisfaction, sort = TRUE)
```

```{r}
# Labling 1 into satisfied, 2 into neutral and 3 into dissatisfied
data <- data |>
  mutate(Service_satisfaction = case_when(
    Service_satisfaction == 1 ~ "Satisfied",
    Service_satisfaction == 2 ~ "Neutral",
    Service_satisfaction == 3 ~ "Dissatisfied"
  ))
data$Service_satisfaction <- as.factor(data$Service_satisfaction)
data |> count(Service_satisfaction, sort = TRUE)
```

```{r}
# Reason for satisfaction
## Converting to factor
data <- data |>
  mutate(across(
    c(Sat_friendly_staff, Sat_service_quality, Sat_short_waiting, Sat_service_available, Sat_others),
    ~ factor(.)
  ))
```

```{r}
## Table for service satisfaction
rm_covsum(
  data = data,
  covs = c('Sat_friendly_staff', 'Sat_service_quality', 'Sat_short_waiting', 'Sat_service_available', 'Sat_others'),
  show.tests = TRUE
)
```


```{r}
# Reasons for dissatisfaction with the vaccination service
## Converting to factor
data <- data |>
  mutate(across(
    c(Dis_staff_not_friendly, Dis_dverse_effect, Dis_long_waiting, Dis_service_close, Dis_no_vaccine),
    ~ factor(.)
  ))
```

```{r}
## Table for reasons for dissatisfaction with the vaccination service
rm_covsum(
  data = data,
  covs = c('Dis_staff_not_friendly', 'Dis_dverse_effect', 'Dis_long_waiting', 'Dis_service_close', 'Dis_no_vaccine'),
  show.tests = TRUE
  )
```


```{r}
# Reasons for not vaccinating the child

## Converting to factor
data <- data |>
  mutate(across(
    c(Reasons_no_no_HF, Reasons_no_dont_no_where, Reasons_no_not_important, Reasons_no_hurt_child, Reasons_no_rumor_danger, Reasons_no_religion_culture, Reasons_no_not_know_benefit, Reasons_no_other), 
    ~ factor(.)
  ))
```

```{r}
## Table for reasons for not vaccinating the child
rm_covsum(
  data = data,
  covs = c('Reasons_no_no_HF', 'Reasons_no_dont_no_where', 'Reasons_no_not_important', 'Reasons_no_hurt_child', 'Reasons_no_rumor_danger', 'Reasons_no_religion_culture', 'Reasons_no_not_know_benefit', 'Reasons_no_other'),
  show.tests = TRUE
)
```

```{r}
# Reasons for under vaccinating the child

## Converting to factor
data <- data |>
  mutate(across(
    c(Reasons_under_far, Reasons_under_bad_schedule, Reasons_under_postponed, Reasons_under_no_staff, Reasons_under_closed, Reasons_under_dont_know_return, Reasons_under_others),
    ~ factor(.)
  ))
```

```{r}
## Table for reasons for under vaccinating the child
rm_covsum(
  data = data,
  covs = c('Reasons_under_far', 'Reasons_under_bad_schedule', 'Reasons_under_postponed', 'Reasons_under_no_staff', 'Reasons_under_closed', 'Reasons_under_dont_know_return', 'Reasons_under_others'),
  show.tests = TRUE
)
```