Turkey rearing

A short project i did for a client.

library setup

library(tidyverse)
library(ggthemes)
library(haven)
library(gtsummary)
my_data<-readxl::read_xlsx("Dissertation_Data.xlsx")
view(my_data) 
names(my_data)
 [1] "...1"                                 
 [2] "Sex"                                  
 [3] "Age"                                  
 [4] "Occupation"                           
 [5] "Constraints"                          
 [6] "...6"                                 
 [7] "Slaughtered for consumption"          
 [8] "Period of selling"                    
 [9] "Reasons for different periods of sale"
[10] "Breed"                                
[11] "Reason  for keeping the breed"        
[12] "Rearing system"                       
[13] "Type of Housing"                      
[14] "Source of feed"                       
[15] "Problems related to feed"             
[16] "Common diseases"                      
[17] "Type of treatment administered"       
[18] "Number of turkeys"                    
[19] "Challenges with veterinary drugs"     
[20] "Mortality causes"

Demographics

Sex

plt_df <- my_data %>%
  group_by(Sex) %>%
  mutate(x.lab=paste0(Sex,
                      "\n",
                      "(n=",
                      n(),
                      ")")) |> 
  group_by(x.lab) |> 
  tally() |> 
  mutate(pct = n/sum(n)*100,
         lbl = sprintf("%.1f%%", pct)) |> 
  mutate(variable="variable") 
  

plt_df %>%
 ggplot() +
  aes(x = fct_reorder(x.lab,pct),y=pct,fill=variable) +
  geom_col(width=0.5) +
  geom_text(aes(label=lbl),position=position_stack(vjust=0.5))+
  labs(title = "", y = "percentage", x ="",fill="") +
  ggthemes::scale_fill_economist() +
  ggthemes::theme_tufte() +
  theme(legend.position = "none",
        axis.text.x = element_text(face="bold", 
                                   color="black",
                                   size=8))
ggsave("Sex.png")

Rearing system

plt_df <- my_data %>%
  group_by(`Rearing system`) %>%
  mutate(x.lab=paste0(`Rearing system`,
                      "\n",
                      "(n=",
                      n(),
                      ")")) |> 
  group_by(x.lab) |> 
  tally() |> 
  mutate(pct = n/sum(n)*100,
         lbl = sprintf("%.1f%%", pct)) |> 
  mutate(variable="variable") 
  

plt_df %>%
 ggplot() +
  aes(x = fct_reorder(x.lab,pct),y=pct,fill=variable) +
  geom_col(width=0.5) +
  geom_text(aes(label=lbl),position=position_stack(vjust=0.5))+
  labs(title = "", y = "percentage", x ="",fill="") +
  ggthemes::scale_fill_economist() +
  ggthemes::theme_tufte() +
  theme(legend.position = "none",
        axis.text.x = element_text(face="bold", 
                                   color="black",
                                   size=8))
ggsave("Rearing_system.png")

Exploration

plt_df <- my_data %>%
  group_by(Constraints) %>%
  mutate(x.lab=paste0(Constraints,
                      "\n",
                      "(n=",
                      n(),
                      ")")) |> 
  group_by(x.lab) |> 
  tally() |> 
  mutate(pct = n/sum(n)*100,
         lbl = sprintf("%.1f%%", pct)) |> 
  mutate(variable="variable") 
plt_df %>%
 ggplot() +
  aes(x = fct_reorder(x.lab,pct),y=pct,fill=variable) +
  geom_col(width=0.5) +
  geom_text(aes(label=lbl),position=position_stack(vjust=0.5))+
  labs(title = "", y = "percentage", x ="",fill="") +
  ggthemes::scale_fill_gdocs() +
  ggthemes::theme_tufte() +
  coord_flip()+
  theme(legend.position = "none",
        plot.title = element_text(face = "bold"),
        axis.text.x = element_text(size = 12),
        axis.ticks.x = element_blank())
ggsave("Constraints.png")

Type of Housing

temp.2 = my_data %>%
  group_by(`Type of Housing`) %>% 
  summarise(count = n()) %>% 
  mutate(per = round(count/sum(count)*100, 2))

ggplot(data = temp.2,
       aes(x = "",
           y = count,
           fill = `Type of Housing`)) +
  geom_col(position = "fill") +
  coord_polar("y",
              start = 10) +
  geom_text(aes(label = paste0(per, "%")),
            position = position_fill(vjust = 0.5),
            size = 4) +
  ggthemes::scale_fill_tableau()+
  ggthemes::theme_tufte()+
  theme(axis.text.x = element_blank(),
        panel.grid = element_blank(),
        axis.ticks = element_blank(),
        legend.position = "right") +
  labs(title = "",
       fill = "Type of Housing",
       x = NULL,
       y = NULL)
ggsave("housing_type.png")

Breed

temp.2 = my_data %>%
  group_by(Breed) %>% 
  summarise(count = n()) %>% 
  mutate(per = round(count/sum(count)*100, 2))

ggplot(data = temp.2,
       aes(x = "",
           y = count,
           fill = Breed)) +
  geom_col(position = "fill") +
  coord_polar("y",
              start = 10) +
  geom_text(aes(label = paste0(per, "%")),
            position = position_fill(vjust = 0.5),
            size = 4) +
  ggthemes::scale_fill_tableau()+
  ggthemes::theme_tufte()+
  theme(axis.text.x = element_blank(),
        panel.grid = element_blank(),
        axis.ticks = element_blank(),
        legend.position = "right") +
  labs(title = "",
       fill = "breed",
       x = NULL,
       y = NULL)
ggsave("breed.png")

Mortality causes

temp.2 = my_data %>%
  group_by(`Mortality causes`) %>% 
  summarise(count = n()) %>% 
  mutate(per = round(count/sum(count)*100, 2))

ggplot(data = temp.2,
       aes(x = "",
           y = count,
           fill = `Mortality causes`)) +
  geom_col(position = "fill") +
  coord_polar("y",
              start = 10) +
  geom_text(aes(label = paste0(per, "%")),
            position = position_fill(vjust = 0.5),
            size = 4) +
  ggthemes::scale_fill_ptol()+
  ggthemes::theme_tufte()+
  theme(axis.text.x = element_blank(),
        panel.grid = element_blank(),
        axis.ticks = element_blank(),
        legend.position = "right") +
  labs(title = "",
       fill = "`Mortality causes`",
       x = NULL,
       y = NULL)
ggsave("mortality.png")
names(my_data)
 [1] "...1"                                 
 [2] "Sex"                                  
 [3] "Age"                                  
 [4] "Occupation"                           
 [5] "Constraints"                          
 [6] "...6"                                 
 [7] "Slaughtered for consumption"          
 [8] "Period of selling"                    
 [9] "Reasons for different periods of sale"
[10] "Breed"                                
[11] "Reason  for keeping the breed"        
[12] "Rearing system"                       
[13] "Type of Housing"                      
[14] "Source of feed"                       
[15] "Problems related to feed"             
[16] "Common diseases"                      
[17] "Type of treatment administered"       
[18] "Number of turkeys"                    
[19] "Challenges with veterinary drugs"     
[20] "Mortality causes"

Type of treatment

plt_df <- my_data %>%
  group_by(`Type of treatment administered`) %>%
  mutate(x.lab=paste0(`Type of treatment administered`,
                      "\n",
                      "(n=",
                      n(),
                      ")")) |> 
  group_by(x.lab) |> 
  tally() |> 
  mutate(pct = n/sum(n)*100,
         lbl = sprintf("%.1f%%", pct)) |> 
  mutate(variable="variable") 
  

plt_df %>%
 ggplot() +
  aes(x = fct_reorder(x.lab,pct),y=pct,fill=variable) +
  geom_col(width=0.5) +
  geom_text(aes(label=lbl),position=position_stack(vjust=0.5))+
  labs(title = "", y = "percentage", x ="",fill="") +
  ggthemes::scale_fill_tableau() +
  ggthemes::theme_tufte() +
  theme(legend.position = "none",
        axis.text.x = element_text(face="bold", 
                                   color="black",
                                   size=8))

Common diseases

my_data |> 
  select(`Common diseases`) |> 
  mutate_if(is.character,as.factor) |> 
  tbl_summary(statistic = list(all_continuous()~"{mean}({sd})",
                               all_categorical()~"{p}%({n}/{N})"))
Characteristic N = 201
Common diseases
    avian influenza 5.0%(1/20)
    foulpox 5.0%(1/20)
    foulpox and infectious coryza 10%(2/20)
    infectious coryza 50%(10/20)
    new castle disease 10%(2/20)
    none 20%(4/20)
1 %(n/N) 
plt_df <- my_data %>%
  group_by(`Problems related to feed`) %>%
  mutate(x.lab=paste0(`Problems related to feed`,
                      "\n",
                      "(n=",
                      n(),
                      ")")) |> 
  group_by(x.lab) |> 
  tally() |> 
  mutate(pct = n/sum(n)*100,
         lbl = sprintf("%.1f%%", pct)) |> 
  mutate(variable="variable") 
  

plt_df %>%
 ggplot() +
  aes(x = fct_reorder(x.lab,pct),y=pct,fill=variable) +
  geom_col(width=0.5) +
  geom_text(aes(label=lbl),position=position_stack(vjust=0.5))+
  labs(title = "", y = "percentage", x ="",fill="") +
  ggthemes::scale_fill_economist() +
  ggthemes::theme_tufte() +
  theme(legend.position = "none",
        axis.text.x = element_text(face="bold", 
                                   color="black",
                                   size=8))
ggsave("problems_related_feed.png")
mod<-aov(`Number of turkeys`~`Common diseases`,data=my_data)
pander::pander(mod)
Analysis of Variance Model
  Df Sum Sq Mean Sq F value Pr(>F)
Common diseases 5 647.4 129.5 2.179 0.1152
Residuals 14 831.6 59.4 NA NA 
mod<-aov(`Number of turkeys`~Constraints,data=my_data)
pander::pander(mod)
Analysis of Variance Model
  Df Sum Sq Mean Sq F value Pr(>F)
Constraints 3 142.9 47.64 0.5705 0.6425
Residuals 16 1336 83.51 NA NA 
my_data %>%
  tbl_summary(include = c(Constraints, 
                          Breed,
                          `Rearing system`,
                          `Mortality causes`,
                          `Rearing system`), 
              by =Constraints,
              type=all_categorical()~"categorical",
              statistic = list(all_continuous()~"{mean}({sd})",
                               all_categorical()~"{p}%({n}/{N})")) %>%
  add_p() %>%
  bold_labels()
Characteristic diseases, N = 81 feed, N = 51 grazing fields, N = 41 predation, N = 31 p-value2
Breed 0.8
    broad breasted bronze 0%(0/8) 20%(1/5) 0%(0/4) 0%(0/3)
    local 88%(7/8) 80%(4/5) 75%(3/4) 100%(3/3)
    mixed 13%(1/8) 0%(0/5) 25%(1/4) 0%(0/3)
Rearing system >0.9
    intensive 50%(4/8) 40%(2/5) 50%(2/4) 33%(1/3)
    semi intensive 50%(4/8) 60%(3/5) 50%(2/4) 67%(2/3)
Mortality causes 0.6
    diseases 25%(2/8) 40%(2/5) 50%(2/4) 0%(0/3)
    none 38%(3/8) 0%(0/5) 25%(1/4) 0%(0/3)
    predation 25%(2/8) 20%(1/5) 0%(0/4) 67%(2/3)
    predation and diseases 13%(1/8) 40%(2/5) 25%(1/4) 33%(1/3)
1 %(n/N)
2 Fisher’s exact test 
my_data %>%
  tbl_summary(include = c(Constraints, 
                          `Problems related to feed`,
                          `Type of Housing`), 
              by =Constraints,
              type=all_categorical()~"categorical",
              statistic = list(all_continuous()~"{mean}({sd})",
                               all_categorical()~"{p}%({n}/{N})")) %>%
  add_p() %>%
  bold_labels()
Characteristic diseases, N = 81 feed, N = 51 grazing fields, N = 41 predation, N = 31 p-value2
Problems related to feed 0.5
    expensive 75%(6/8) 60%(3/5) 50%(2/4) 67%(2/3)
    inadequate knowledge on feed 25%(2/8) 40%(2/5) 50%(2/4) 0%(0/3)
    no problem 0%(0/8) 0%(0/5) 0%(0/4) 33%(1/3)
Type of Housing 0.2
    mash wire 38%(3/8) 40%(2/5) 0%(0/4) 33%(1/3)
    mash wire and asbestos 0%(0/8) 0%(0/5) 0%(0/4) 33%(1/3)
    mash wire and bricks 38%(3/8) 60%(3/5) 100%(4/4) 0%(0/3)
    mash wire and old metal 13%(1/8) 0%(0/5) 0%(0/4) 0%(0/3)
    mash wire and wood 13%(1/8) 0%(0/5) 0%(0/4) 33%(1/3)
1 %(n/N)
2 Fisher’s exact test 
my_data |> 
  select(Sex,Age,Occupation) |> 
  mutate_if(is.character,as.factor) |> 
  tbl_summary(statistic = list(all_continuous()~"{mean}({sd})",
                               all_categorical()~"{p}%({n}/{N})"))
Characteristic N = 201
Sex
    female 70%(14/20)
    male 30%(6/20)
Age 36(10)
Occupation
    family member 25%(5/20)
    owner 75%(15/20)
1 %(n/N); Mean(SD) 
my_data |> 
  select(`Reason  for keeping the breed`) |> 
  mutate_if(is.character,as.factor) |> 
  tbl_summary(statistic = list(all_continuous()~"{mean}({sd})",
                               all_categorical()~"{p}%({n}/{N})"))
Characteristic N = 201
Reason for keeping the breed
    adaptability 10%(2/20)
    affordable 30%(6/20)
    availability 25%(5/20)
    fast growth 5.0%(1/20)
    genetic variation 5.0%(1/20)
    inadequate knowledge 5.0%(1/20)
    resistant to diseases 10%(2/20)
    versatile diet 10%(2/20)
1 %(n/N) 
mod<-t.test(`Number of turkeys`~`Rearing system`,data=my_data)
pander::pander(mod)
Welch Two Sample t-test: Number of turkeys by Rearing system (continued below)
Test statistic df P value Alternative hypothesis
-2.662 17.89 0.01593 * two.sided
mean in group intensive mean in group semi intensive
8.556 17.55 
names(data)
NULL
NicotineL~gender+education+Oscore+Cscore+Impulsive
NicotineL ~ gender + education + Oscore + Cscore + Impulsive