#DATA CLEANING

#Load libraries ---- 

library(tidyverse)
library(dplyr)
library(readxl)
library(janitor)
library(here)

#Read in data ---- 

ipa <- read_excel("Krisengefüge-Einfluss_der_Corona_Pandemie_auf_die_darstellenden_Künste_in_den_D-A-CH_Ländern-4.xlsx")

#Tidy data ---- 

rename_ipa <- rename(ipa,
                     role = A001, 
                     age = QF01_01,
                     sex = QF02,
                     
                     nationality_a = CODE_QF13_01,
                     nationality_b = CODE_b_QF13_01,
                     location =  QF12,
                     
                     degree = QF04,
                     employment = QF07,
                     
                     nr_prof = QF05,
                     artist = QF05_01,
                     designer = QF05_02,
                     technician = QF05_03,
                     producer = QF05_04,
                     dramaturge =  QF05_05,
                     curator = QF05_06,
                     art_manager = QF05_07,
                     other_prof = QF05_08,
                     other_prof_name = CODE_QF05_08a,
                     
                     nr_genre = QF06,
                     dance = QF06_01,
                     theater = QF06_02,
                     performance = QF06_03,
                     music_theater = QF06_04,
                     child_youth = QF06_05,
                     circus = QF06_06,
                     interdis = QF06_07,
                     comercial = QF06_10,
                     ad_pr = QF06_09,
                     other_genre = QF06_08, 
                     other_genre_name = CODE_QF06_08a,
                     
                     nr_soci_sec = QF09,
                     publ_soci_sec = QF09_01,
                     artist_soci_sec = QF09_03,
                     artist_soci_sec_name = CODE_QF09_03a,
                     no_soci_sec = QF09_02,
                     na_social_sec = QF09_04,
                     
                     member = QF10,
                     member_eaipa = CODE_QF10_03,
                     member_other = CODE_b_QF10_03,
                     
                     country = CODE_QF11_01,
                     province = CODE_QF11_02,
                     city = CODE_QF11_03,
                     residence_b = CODE_b_location,
                     
                     nr_clients19 = QF14_01,
                     nr_clients20 = QF14_02,
                     
                     need_f_com = TO01_01,
                     need_f_collab = TO01_02,
                     pol_alliance = TO01_03,
                     ident_w_ipa = TO01_04,
                     incr_info = TO01_05,
                     build_struc = TO01_06,
                     joint_venture = TO01_07,
                     
                     currency = CODE_SO06_01,
                     wage_nominal = SO01,
                     wage = SO01_01,
                     wage_eur = CODE_SO01_01,
                     income19 = SO02_01,
                     income19_eur = CODE_SO02_01,
                     income20 = SO02_03,
                     income20_eur = CODE_SO02_03,
                     
                     unempl_insur = SO03_01,
                     unempl_corona = SO03_02,
                     gov_sup_corona = SO03_03,
                     empl_ouside_ipa = SO09_01,
                     income_plenty = SO10_01,
                     income_adequacy = SO10_03,
                     job_stabi = S103_01,
                     fin_secu = S103_02,
                     pension_secu = S103_03,
                     unpaid_work = S103_04,
                     pos_toward_overwork = S103_05,
                     ident_w_prof = S103_06,
                     meaningfullness_o_work = S103_07,
                     
                     rep_through_advo = AA01_01,
                     particip_i_advo = AA01_02,
                     policy_feedback_pos = AA01_03,
                     policy_feedback_consid = AA01_04,
                     dialog_gov_mun = AA01_05,
                     dialog_gov_nat = AA01_06,
                     policy_com_gov = AA01_07,
                     demo_asso = AA01_09,
                     coop_gov = AA01_10,
                     policy_inclu = AA01_11,
                     
                     nr_orgatype = AY01,
                     company = AY01_01,
                     product_house = AY01_02,
                     rehear_space = AY01_03,
                     venue_w_ensam = AY01_04,
                     venue_wo_ensam = AY01_05,
                     festival = AY01_06,
                     resid_space = AY01_08,
                     orgatype_other = AY01_07,
                     orgatype_other_name = CODE_AY01_07a,
                     
                     npo_stat = AY02_01,
                     fin_indep = AY02_02,
                     legal_stat = AY03,
                     legal_stat_other_name = AY03_06,
                     
                     member_org = AY04,
                     member_eaipa_org = CODE_AY04_03,
                     member_other_org = CODE_b_AY04_03,
                     
                     nr_genre_org = AY05,
                     dance_org = AY05_01,
                     theater_org = AY05_02,
                     performance_org = AY05_03,
                     music_theater_org = AY05_04,
                     child_youth_org = AY05_05,
                     circus_org = AY05_06,
                     interdis_org = AY05_07,
                     other_genre_org = AY05_08,
                     other_genre_org_name = CODE_AY05_08a,
                     
                     country_org = CODE_AY06_01,
                     province_org = CODE_AY06_02,
                     city_org = CODE_AY06_03,
                     location_org = AY07,
                     
                     need_f_com_org = OF01_01,
                     need_f_collab_org = OF01_02,
                     pol_alliance_org =OF01_03,
                     ident_w_ipa_org = OF01_04,
                     incr_info_org = OF01_05,
                     build_struc_org = OF01_06,
                     joint_venture_org = OF01_07,
                     
                     full_time19 = S201_02,
                     part_time19 = S201_03,
                     temporary19 = S201_04,
                     fee19 = S201_05,
                     special_order19 = S201_01,
                     volun_intern19 = S201_06,
                     no_cot19 = S201_07,
                     full_time20 = S203_02,
                     part_time20 = S203_03,
                     temporary20 = S203_04,
                     fee20 = S203_05,
                     special_order20 = S203_01,
                     volun_intern20 = S203_06,
                     no_cot20 = S203_07,
                     
                     rep_through_advo_org = A101_01,
                     particip_i_advo_org = A101_02,
                     policy_feedback_pos_org = A101_03,
                     policy_feedback_consid_org = A101_04,
                     dialog_gov_mun_org = A101_05,
                     dialog_gov_nat_org = A101_06,
                     policy_com_gov_org = A101_07,
                     demo_asso_org = A101_08,
                     coop_gov_org = A101_09,
                     policy_inclu_org = A101_10,
                     
                     prem19 = OO01_01,
                     prod19 = OO01_02,
                     shows_avre19 = OO01_06,
                     shows19 = OO01_07,
                     guestperf19 = OO01_08,
                     spectators19 = OO01_13,
                     
                     prem20 = OO06_01,
                     prod20 = OO06_02,
                     shows_avr20 = OO06_03,
                     shows20 = OO06_04,
                     guestperf20 = OO06_05,
                     spectators20 = OO06_07,
                     
                     own_venue = OO02,
                     audcapa = OO02_01,
                     own_rehear_space = OO03,
                     sice = OO03_01,
                     
                     planing_secu_org = OO04,
                     pub_fund_org = OO07,
                     
                     comment = F102_01)

clean_ipa <- clean_names(rename_ipa)
names(clean_ipa)

# Sumarize to check new names 
summary(clean_ipa)
glimpse(clean_ipa)

# Compute total variables and change all -9 to NA ---- 

clean_ipa <- clean_ipa %>% mutate(location = coalesce(location,location_org),
                      member = coalesce(member, member_org),
                      member_eaipa = coalesce(member_eaipa, member_eaipa_org),
                      member_other = coalesce(member_other, member_other_org),
                      country = coalesce(country, country_org),
                      state_province = coalesce(province, province_org),
                      city = coalesce(city, city_org),
                      need_f_com = coalesce(need_f_com, need_f_com_org),
                      need_f_collab = coalesce(need_f_collab, need_f_collab_org),
                      pol_alliance = coalesce(pol_alliance, pol_alliance_org),
                      ident_w_ipa = coalesce(ident_w_ipa, ident_w_ipa_org),
                      incr_info = coalesce(incr_info, incr_info_org),
                      build_struc = coalesce(build_struc, build_struc_org),
                      joint_venture = coalesce(joint_venture, joint_venture_org),
                      rep_through_advo = coalesce(rep_through_advo, rep_through_advo_org),
                      particip_i_advo = coalesce(particip_i_advo, particip_i_advo_org),
                      policy_feedback_pos = coalesce(policy_feedback_pos, policy_feedback_pos_org),
                      policy_feedback_consid = coalesce(policy_feedback_consid, policy_feedback_consid_org),
                      dialog_gov_mun = coalesce(dialog_gov_mun, dialog_gov_mun_org),
                      dialog_gov_nat = coalesce(dialog_gov_nat, dialog_gov_nat_org),
                      policy_com_gov = coalesce(policy_com_gov, policy_com_gov_org),
                      demo_asso = coalesce(demo_asso, demo_asso_org),
                      coop_gov = coalesce(coop_gov, coop_gov_org),
                      policy_inclu = coalesce(policy_inclu, policy_inclu_org),
                                  
                      nr_genre = coalesce(nr_genre, nr_genre_org),
                      dance = coalesce(dance, dance_org),
                      theater = coalesce(theater, theater_org),
                      performance = coalesce(performance, performance_org),
                      music_theater = coalesce(music_theater, music_theater_org),
                      child_youth = coalesce(child_youth, child_youth_org),
                      circus = coalesce(circus, circus_org),
                      interdis = coalesce(interdis, interdis_org),
                                
                      other_genre = coalesce(other_genre, other_genre_org),
                      other_genre_name = coalesce(other_genre_name, other_genre_org_name))

# tralslate all -9 to na and all -1 to 0
clean_ipa[clean_ipa == -9] <- NA
clean_ipa[clean_ipa == -1] <- NA

# Combine variables ---- 

clean_ipa$eaipa_country <- clean_ipa$country
clean_ipa$eaipa_country[clean_ipa$eaipa_country %in% c("United Kingdom",
                      "The Netherlands",
                      "Spain",
                      "Slovakia",
                      "Portugal",
                      "Norway",
                      "Croatia",
                      "Belgium",
                      "ù")] <- "Other"

clean_ipa$eaipa_country[clean_ipa$eaipa_country %in% c("Switzerland",
                      "Liechtenstein")] <- "Switzerland"

clean_ipa$eu_regions[clean_ipa$country %in% c("Austria",
                      "Germany",
                      "Liechtenstein",
                      "Switzerland",
                      "The Netherlands",
                      "Belgium")] <- "Western Europe"

clean_ipa$eu_regions[clean_ipa$country %in% c("Bulgaria",
                      "Croatia",
                      "Czech Republic",
                      "Hungary",
                      "Romania",
                      "Slovakia")] <- "Eastern Europe"

clean_ipa$eu_regions[clean_ipa$country %in% c("Spain",
                      "Italy",
                      "Slovenia")] <- "Southern Europe"

clean_ipa$eu_regions[clean_ipa$country %in% c("Finland",
                      "Iceland",
                      "Norway",
                      "Sweden",
                      "United Kingdom")] <- "Northern Europe"

# Create Factos ---- 

sum_ipa <- clean_ipa %>% 
  mutate(
    sex = factor(sex,
                      levels = c("1", "2", "3"),
                      labels = c("Female", "Male", "other")),
    role = factor(role,
                      levels = c("1", "2"),
                      labels = c("organization", "individual")),
    location = factor(location,
                      levels = c("1", "2", "3", "4"),
                      labels = c("large_city", "city", "town", "rural")),
    degree = factor(degree,ordered = TRUE),
    employment = factor(employment, levels = c("1", "2", "3", "4"),
                      labels = c("freelancer", "part-time employee", "full-time employee", "hybrid (freelance/empl)")),
    publ_soci_sec = factor(publ_soci_sec, ordered = TRUE),
    artist_soci_sec = factor(artist_soci_sec, ordered = TRUE),
    no_soci_sec = factor(no_soci_sec, ordered = TRUE),
    na_social_sec = factor(na_social_sec, ordered = TRUE),
    member = factor(member,
                      levels = c("2", "3"),
                      labels = c("no", "yes")),
    
    nr_clients19 = factor(nr_clients19,ordered = TRUE),
    nr_clients20 = factor(nr_clients20,ordered = TRUE),
    
    unempl_insur = factor(unempl_insur, ordered = TRUE),
    unempl_corona = factor(unempl_corona, ordered = TRUE),
    gov_sup_corona = factor(gov_sup_corona, ordered = TRUE),
    empl_ouside_ipa = factor(empl_ouside_ipa,ordered = TRUE),
    npo_stat = factor(npo_stat,
                      levels = c("1", "2"),
                      labels = c("No", "Yes")), 
    fin_indep = factor(fin_indep,
                      levels = c("1", "2"),
                      labels = c("No", "Yes")), 
    legal_stat = factor(legal_stat,
                      levels = c("1", "2", "3", "4", "5", "6"),
                      labels = c("Private company",
                      "Limited liability company",
                      "Association",
                      "Cooperative",
                      "Foundation",
                      "Other")), 
    location_org = factor(location_org,
                      levels = c("1", "2", "3", "4"),
                      labels = c("large_city", "city", "town", "rural")),
    member_org = factor(member_org,
                      levels = c("2", "3"),
                      labels = c("no", "yes")),
    own_venue = factor(own_venue,
                      levels = c("1", "2"),
                      labels = c("Yes", "No")), 
    own_rehear_space = factor(own_rehear_space,
                      levels = c("1", "2"),
                      labels = c("Yes", "No")), 
    planing_secu_org = factor(planing_secu_org,ordered = TRUE, 
                      levels = c("7", "1", "2", "3", "4", "5", "6")),
    pub_fund_org = factor(pub_fund_org,ordered = TRUE, 
                      levels = c("5", "4", "3", "2", "1")),
    
    need_f_com = factor(need_f_com,ordered = TRUE),
    need_f_collab = factor(need_f_collab,ordered = TRUE),
    pol_alliance = factor(pol_alliance,ordered = TRUE),
    ident_w_ipa = factor(ident_w_ipa,ordered = TRUE),
    incr_info = factor(incr_info,ordered = TRUE),
    build_struc = factor(build_struc,ordered = TRUE),
    joint_venture = factor(joint_venture,ordered = TRUE),
    
    income_plenty = factor(income_plenty,ordered = TRUE),
    income_adequacy = factor(income_adequacy, ordered = TRUE),
    job_stabi = factor(job_stabi,ordered = TRUE),
    fin_secu = factor(fin_secu,ordered = TRUE),
    pension_secu = factor(pension_secu,ordered = TRUE),
    
    unpaid_work = factor(unpaid_work,ordered = TRUE),
    pos_toward_overwork = factor(pos_toward_overwork,ordered = TRUE),
    ident_w_prof = factor(ident_w_prof,ordered = TRUE),
    meaningfullness_o_work = factor(meaningfullness_o_work,ordered = TRUE),
    
    rep_through_advo = factor(rep_through_advo,ordered = TRUE),
    particip_i_advo = factor(particip_i_advo,ordered = TRUE),
    policy_feedback_pos = factor(policy_feedback_pos,ordered = TRUE),
    policy_feedback_consid = factor(policy_feedback_consid,ordered = TRUE),
    dialog_gov_mun = factor(dialog_gov_mun,ordered = TRUE),
    dialog_gov_nat = factor(dialog_gov_nat,ordered = TRUE),
    policy_com_gov = factor(policy_com_gov,ordered = TRUE),
    demo_asso = factor(demo_asso,ordered = TRUE),
    coop_gov = factor(coop_gov,ordered = TRUE),
    policy_inclu = factor(policy_inclu,ordered = TRUE),
    
    need_f_com_org = factor(need_f_com_org,ordered = TRUE),
    need_f_collab_org = factor(need_f_collab_org,ordered = TRUE),
    pol_alliance_org = factor(pol_alliance_org,ordered = TRUE),
    ident_w_ipa_org   = factor(ident_w_ipa_org,ordered = TRUE),
    incr_info_org = factor(incr_info_org,ordered = TRUE),
    build_struc_org = factor(build_struc_org,ordered = TRUE),
    joint_venture_org  = factor(joint_venture_org,ordered = TRUE),
    
    full_time19 = factor(full_time19,ordered = TRUE),
    part_time19 = factor(part_time19,ordered = TRUE),
    temporary19 = factor(temporary19,ordered = TRUE),
    fee19 = factor(fee19,ordered = TRUE),
    special_order19 = factor(special_order19,ordered = TRUE),
    volun_intern19 = factor(volun_intern19,ordered = TRUE),
    no_cot19 = factor(no_cot19,ordered = TRUE),
    full_time20 = factor(full_time20,ordered = TRUE),
    part_time20 = factor(part_time20,ordered = TRUE),
    temporary20 = factor(temporary20,ordered = TRUE),
    fee20 = factor(fee20,ordered = TRUE),
    special_order20 = factor(special_order20,ordered = TRUE),
    volun_intern20 = factor(volun_intern20,ordered = TRUE),
    no_cot20 = factor(no_cot20,ordered = TRUE),
    rep_through_advo_org = factor(rep_through_advo_org,ordered = TRUE),
    particip_i_advo_org = factor(particip_i_advo_org,ordered = TRUE),
    policy_feedback_pos_org = factor(policy_feedback_pos_org,ordered = TRUE),
    policy_feedback_consid_org = factor(policy_feedback_consid_org,ordered = TRUE),
    dialog_gov_mun_org = factor(dialog_gov_mun_org,ordered = TRUE),
    dialog_gov_nat_org = factor(dialog_gov_nat_org,ordered = TRUE),
    policy_com_gov_org = factor(policy_com_gov_org,ordered = TRUE),
    demo_asso_org = factor(demo_asso_org,ordered = TRUE),
    coop_gov_org = factor(coop_gov_org,ordered = TRUE),
    policy_inclu_org = factor(policy_inclu_org,ordered = TRUE))

clean_ipa$eaipa_country <- factor(clean_ipa$eaipa_country)
clean_ipa$country <- factor(clean_ipa$country)
str(sum_ipa)

# separate in usefull dataframes ---- 

ipa_ind <- subset(sum_ipa, role == "individual")
ipa_ind <- ipa_ind[, -90:-201] 
ipa_org <- subset(sum_ipa, role == "organization")
ipa_org <- ipa_org[,-8:-89, -171:-201] 
ipa_common <- sum_ipa[ ,c(1,4:7,12,25:36,43:45,46:48,52:58,80:89,176:182,190,191,193:201, 202, 203)]

# Check plausibility and exclude implausible ---- 

sapply(sapply(ipa_ind, unique), sort, na.last=TRUE)
colSums(is.na(ipa_common))
colSums(is.na(ipa_org))
colSums(is.na(ipa_ind))

summary(ipa_common)
# Inconsistencies are:
  # country_org UNPLAUSIBLE: ù
  # state_province_org UNPLAUSIBLE: ù
  # city_org UNPLAUSIBLE: ù
  # audcapa_venue UNPLAUSIBLE: 1000000
  # in wage there are too big numbers / in icome20 too, 
  # premieres19 CHECK: 123 133

sum_ipa_check <- subset(sum_ipa, select = c(case, country_org, city_org, province_org, audcapa))
# Inconsistencies are:
  # nr 1239 gives only repetitive answers and the ù instead of text fields. 
  # nr 2983 has given improbable capacity
  # nr 466 Income 20 is implausibly high

# Exclude Inconsistencies 
sum_ipa <- sum_ipa[-c(232, 596, 621), ]

# Apply to all subsets
ipa_ind <- subset(sum_ipa, role == "individual")
ipa_ind <- ipa_ind[, -90:-201] 
ipa_org <- subset(sum_ipa, role == "organization")
ipa_org <- ipa_org[,-8:-89, -171:-201] 
ipa_common <- sum_ipa[ ,c(1,4:7,12,25:36,43:45,46:48,52:58,80:89,176:182,190,191,193:201, 202, 203, 204)]

ipa_ind_check <- subset(ipa_ind, select = c(case, artist, designer, technician, producer, dramaturge, curator, art_manager, country, employment, wage_eur, income19_eur, income20_eur))
# Inconsistencies are:
  #nr 287 / 23 wage 13537 EUR implausible, not reflected in income
  #nr 791 / 110 wage 32600 EUR implausible, not reflected in income
  #nr 915 / 131 wage 9000 EUR implausible, not reflected in income
  #no. 2923 / 373 wages 6134 EUR implausible, not reflected in income
  #nr 3185 / 421 wage 12000 EUR implausible, not reflected in income
  #nr 3972 / 509 wage 19560 EUR implausible, not reflected in income

# Exclude Inconsistencies 
ipa_wage <- ipa_ind[-c(23, 110, 131, 372, 321, 509), ]

ipa_org_check <- subset(ipa_org, select = c(case, 
                                            prem19, prod19, shows_avre19, shows19, guestperf19, spectators19, 
                                            prem20, prod20, shows_avr20, shows20, guestperf20, spectators20,
                                            full_time19, part_time19, temporary19, fee19, special_order19, volun_intern19, no_cot19,
                                            full_time20, part_time20, temporary20, fee20, special_order20, volun_intern20, no_cot20))
# Inconsistencies are:
  # 1874, 133 premieres. Implausible numbers --> deleted
  # 2338, 123 premieres. Implausible numbers --> deleted

# Exclude Inconsistencies 
ipa_org <- ipa_org[-c(130, 170), ]

# Compute long dataframes ----
ipa_profession <- ipa_ind %>% pivot_longer(names_to = "profession", values_to = "yes_no", artist:other_prof)
ipa_profession <- filter(ipa_profession, yes_no == 2)

ipa_orgatype <- ipa_org %>% pivot_longer(names_to = "orgatype", values_to = "orgvalue", company:orgatype_other)
ipa_orgatype <- filter(ipa_orgatype, orgvalue == 2)

ipa_genre_all <- ipa_common %>% pivot_longer(names_to = "genre_all", values_to = "yes_no", dance:other_genre)
ipa_genre_all <- filter(ipa_genre_all, yes_no == 2)

ipa_genre_all$genre_all[ipa_genre_all$genre_all %in% c("ad_pr",
                                                       "comercial")] <- "comercial"

# Write new file that can be used for further analysis ----
write_excel_csv(sum_ipa, here("Data", "Data_out", "sum_ipa.csv"))
write_excel_csv(ipa_ind, here("Data", "Data_out", "ipa_ind.csv"))
write_excel_csv(ipa_org, here("Data", "Data_out", "ipa_org.csv"))
write_excel_csv(ipa_common, here("Data", "Data_out", "ipa_common.csv"))
write_excel_csv(ipa_wage, here("Data", "Data_out", "_money", "ipa_wage.csv"))
write_excel_csv(ipa_profession, here("Data", "Data_out", "_prof_org", "ipa_profession.csv"))
write_excel_csv(ipa_orgatype, here("Data", "Data_out", "_prof_org", "ipa_orgatype.csv"))
write_excel_csv(ipa_genre_all, here("Data", "Data_out", "_genre", "ipa_genre_all.csv"))

DATA INSPECTION 

library(tidyverse)
library(dplyr)
library(psych)
library(GGally)

# Read in Data ----

sum_ipa <- read.csv("Data/Data_out/sum_ipa.csv")
ipa_ind <- read.csv("Data/Data_out/ipa_ind.csv")
ipa_org <- read.csv("Data/Data_out/ipa_org.csv")
ipa_common <- read.csv("Data/Data_out/ipa_common.csv")
ipa_profession <- read_csv("Data/Data_out/_prof_org/ipa_profesion.csv")
ipa_genre_all <- read_csv("Data/Data_out/_genre/ipa_genre_all.csv")
ipa_orgatype <- read_csv("Data/Data_out/_prof_org/ipa_orgatype.csv")
ipa_wage <- read_csv("Data/Data_out/_money/ipa_wage.csv")
clean_ipa <- read_csv("Data/Data_out/clean_ipa.csv")

# inspect individual variables / check normal distribution ----

data <- as.Date.factor(sum_ipa$started)
barplot(table(data), xlab = "Date of submission", ylab = "Number of cases")

hist(clean_ipa$age)
hist(clean_ipa$sex)
hist(clean_ipa$location)

hist(clean_ipa$employment)
hist(clean_ipa$nr_prof)
hist(clean_ipa$nr_genre)
hist(clean_ipa$nr_clients19)
hist(clean_ipa$nr_clients20)

hist(ipa_ind$need_f_com)
hist(ipa_org$need_f_com)
hist(ipa_common$need_f_com)

hist(ipa_ind$need_f_collab)
hist(ipa_org$need_f_collab)
hist(ipa_common$need_f_collab)

hist(ipa_ind$pol_alliance)
hist(ipa_org$pol_alliance)
hist(ipa_common$pol_alliance)

hist(ipa_ind$incr_info)
hist(ipa_org$incr_info)
hist(ipa_common$incr_info)

hist(ipa_ind$build_struc)
hist(ipa_org$build_struc)
hist(ipa_common$build_struc)

hist(ipa_ind$joint_venture)
hist(ipa_org$joint_venture)
hist(ipa_common$joint_venture)

hist(clean_ipa$wage_eur)
hist(clean_ipa$income19_eur)
hist(clean_ipa$income20_eur)

hist(ipa_ind$unempl_insur)
hist(ipa_ind$unempl_corona)
hist(ipa_ind$gov_sup_corona)

hist(ipa_ind$income_plenty)
hist(ipa_ind$income_adequacy)
hist(ipa_ind$job_stabi)
hist(ipa_ind$fin_secu)
hist(ipa_ind$pension_secu)

hist(ipa_ind$unpaid_work)
hist(ipa_ind$pos_toward_overwork)

hist(ipa_ind$ident_w_prof)
hist(ipa_ind$meaningfullness_o_work)

hist(ipa_ind$rep_through_advo)
hist(ipa_org$rep_through_advo)
hist(ipa_common$rep_through_advo)

hist(ipa_ind$particip_i_advo)
hist(ipa_org$particip_i_advo)
hist(ipa_common$particip_i_advo)

hist(ipa_ind$policy_feedback_pos)
hist(ipa_org$policy_feedback_pos)
hist(ipa_common$policy_feedback_pos)

hist(ipa_ind$policy_feedback_consid)
hist(ipa_org$policy_feedback_consid)
hist(ipa_common$policy_feedback_consid)

hist(ipa_ind$dialog_gov_mun)
hist(ipa_org$dialog_gov_mun)
hist(ipa_common$dialog_gov_mun)

hist(ipa_ind$dialog_gov_nat)
hist(ipa_org$dialog_gov_nat)
hist(ipa_common$dialog_gov_nat)

hist(ipa_ind$policy_com_gov)
hist(ipa_org$policy_com_gov)
hist(ipa_common$policy_com_gov)

hist(ipa_ind$demo_asso)
hist(ipa_org$demo_asso)
hist(ipa_common$demo_asso)

hist(ipa_ind$coop_gov)
hist(ipa_org$coop_gov)
hist(ipa_common$coop_gov)

hist(ipa_ind$policy_inclu)
hist(ipa_org$policy_inclu)
hist(ipa_common$policy_inclu)

hist(ipa_org$npo_stat)
names(sum_ipa)

ipa_profession <- filter(ipa_profession, yes_no == 2)

ipa_ind <- filter(ipa_ind, income19_eur < 100000)
ipa_ind <- filter(ipa_ind, income20_eur < 100000)

boxplot(ipa_ind$wage_eur~ipa_ind$sex)
boxplot(ipa_ind$income19_eur~ipa_ind$sex)
boxplot(ipa_ind$income20_eur~ipa_ind$sex)

boxplot(ipa_ind$wage_eur~ipa_ind$degree)
boxplot(ipa_ind$income19_eur~ipa_ind$degree)
boxplot(ipa_ind$income20_eur~ipa_ind$degree)

boxplot(ipa_ind$wage_eur~ipa_ind$country)
boxplot(ipa_ind$income19_eur~ipa_ind$country)
boxplot(ipa_ind$income20_eur~ipa_ind$country)

#Scatterplot Matrix to see if indices correlate

pairs(~ need_f_com + need_f_collab + pol_alliance + ident_w_ipa + incr_info + build_struc + joint_venture, data = ipa_common,
      lower.panel = panel.smooth)

sum_ipa %>% ggpairs(columns = c("need_f_com",
                                  "need_f_collab",
                                  "pol_alliance",
                                  "ident_w_ipa",
                                  "incr_info",
                                  "build_struc",
                                  "joint_venture"))

ipa_ind %>% ggpairs(columns = c("publ_soci_sec",
                                 "no_soci_sec",
                                 "unempl_insur"))
                                 
ipa_ind %>% ggpairs(columns = c("income_plenty",    
                                "income_adequacy",  
                                "fin_secu", 
                                "job_stabi",
                                "pension_secu",
                                "unpaid_work"))

ipa_common %>% ggpairs(columns = c("structuration_scale",    
                                "governance_scale"
                                ))

# SAMPEL DESCRIBTION ALL  ---------------

tabyl(ipa_common, role, )
tabyl(ipa_ind, role, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Survey participant distribution by country
tabyl(ipa_common, country, )
tabyl(ipa_ind, country, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Survey participant distribution by ind org
tabyl(sum_ipa, role, country)
tabyl(sum_ipa , role, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Distribution genre 
ipa_genre_all <- ipa_common %>% pivot_longer(names_to = "genre_all", values_to = "yes_no", dance:other_genre)
ipa_genre_all <- filter(ipa_genre_all, yes_no == 2)

tabyl(ipa_genre_all, genre_all, )
tabyl(ipa_genre_all, genre_all, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Number genre 
tabyl(ipa_common, nr_genre, ) 
tabyl(ipa_common, nr_genre, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Location
tabyl(ipa_common,location, )
tabyl(ipa_common,location, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Member
tabyl(ipa_common,member, )
tabyl(ipa_common,member, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Member EAIPA
tabyl(ipa_common,member_eaipa, country)
tabyl(ipa_common,member_eaipa, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# SAMPEL DESCRIBTION INDIVIDUALS  ----

# Gender
tabyl(ipa_ind,sex, )
tabyl(ipa_ind,sex, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Age 
mean(ipa_ind$age)
describe(ipa_ind$age)
describeBy(ipa_ind$age, ipa_ind$country) 

# Nationality
tabyl(ipa_ind, nationality_a, )
tabyl(ipa_ind, nationality_a, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Education
tabyl(ipa_ind,degree, location) 
tabyl(ipa_ind, degree) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

describe.by(ipa_ind$degree, ipa_ind$country)

# Profession
ipa_profession$country <- ipa_profession$country
ipa_profession$country[ipa_profession$country %in% c("United Kingdom",
                                                       "The Netherlands",
                                                       "Spain",
                                                       "Slovakia",
                                                       "Portugal",
                                                       "Norway",
                                                       "Croatia",
                                                       "Belgium",
                                                       "ù")] <- "Other"
ipa_profession$country[ipa_profession$country %in% c("Switzerland",
                                                       "Liechtenstein")] <- "Switzerland"

ipa_profession <- ipa_profession %>% pivot_longer(names_to = "profession", values_to = "yes_no", artist:other_prof)
ipa_profession <- filter(ipa_profession, yes_no == 2)

# Distribution Profession
tabyl(ipa_profession, profession, country) 
tabyl(ipa_profession, profession, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Number Profession
tabyl(ipa_ind,nr_prof, ) 
tabyl(ipa_ind,nr_prof, )  %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# SAMPEL DESCRIBTION ORGANIZATIONS ----------------------

#Orgatype
ipa_orgatype <- ipa_org %>% pivot_longer(names_to = "orgatype", values_to = "orgvalue", company:orgatype_other)
ipa_orgatype <- filter(ipa_orgatype, orgvalue == 2)

tabyl(ipa_org, nr_orgatype, country_org)
tabyl(ipa_org, nr_orgatype, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Distribution orgatype
tabyl(ipa_orgatype, orgatype, )
tabyl(ipa_orgatype, orgatype, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Legal status
tabyl(ipa_org, legal_stat, )
tabyl(ipa_org, legal_stat, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# NPO status
tabyl(ipa_org, npo_stat, )
tabyl(ipa_org, npo_stat, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# SAMPEL DISCRIBTION PER COUNTRY - ----------------

# Distribution genre 
ipa_genre_all <- ipa_common %>% pivot_longer(names_to = "genre_all", values_to = "yes_no", dance:other_genre)
ipa_genre_all <- filter(ipa_genre_all, yes_no == 2)

tabyl(ipa_genre_all, genre_all, country)
tabyl(ipa_genre_all, genre_all, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Number genre 
tabyl(ipa_common, nr_genre, country) 
tabyl(ipa_common, nr_genre, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Location
tabyl(ipa_common,location, country)
tabyl(ipa_common,location, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Gender
tabyl(ipa_ind, sex, country)
tabyl(ipa_ind, sex, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Age 
tabyl(ipa_ind, age, country)
tabyl(ipa_ind, age, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Education
tabyl(ipa_ind,degree, country) 
tabyl(ipa_ind,degree, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Profession
ipa_profession <- ipa_ind %>% pivot_longer(names_to = "profession", values_to = "yes_no", artist:other_prof)
ipa_profession <- filter(ipa_profession, yes_no == 2)

# Distribution Professions 
tabyl(ipa_profession, profession, country) 
tabyl(ipa_profession, profession, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Number professions
tabyl(ipa_ind,nr_prof, country) 
tabyl(ipa_ind,nr_prof, country)  %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

#Orgatype
ipa_orgatype <- ipa_org %>% pivot_longer(names_to = "orgatype", values_to = "orgvalue", company:orgatype_other)
ipa_orgatype <- filter(ipa_orgatype, orgvalue == 2)

tabyl(ipa_org, nr_orgatype, country)
tabyl(ipa_org, nr_orgatype, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Distribution Orgatype
tabyl(ipa_orgatype, orgatype, country)
tabyl(ipa_orgatype, orgatype, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Legal status
tabyl(ipa_org, legal_stat, country)
tabyl(ipa_org, legal_stat, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# NPO status
tabyl(ipa_org, npo_stat, country)
tabyl(ipa_org, npo_stat, country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# 1. ISTITUTIONALIZATION

#Load libraries ---- 

library(psych)
library(CGPfunctions)
library(DescTools)
library(readr)
library(GGally)
library(mice)

#Read in data ---- 

ipa_common <- read_csv("Data/Data_out/ipa_common.csv")
ipa_ind <- read_csv("Data/Data_out/ipa_ind.csv")
ipa_org <- read_csv("Data/Data_out/ipa_org.csv")
ipa_genre_all <- read_csv("Data/Data_out/_genre/ipa_genre_all.csv")
ipa_profession <- read_csv("Data/Data_out/_prof_org/ipa_profession.csv")
ipa_orgatype <- read_csv("Data/Data_out/_prof_org/ipa_orgatype.csv")

# Genre ----
# Genre group allocation

ipa_common$nr_genreSum[ipa_common$nr_genre %in% c("1",
                                               "0")] <- "1 Genre"
ipa_common$nr_genreSum[ipa_common$nr_genre %in% c("2")] <- "2 Genre"
ipa_common$nr_genreSum[ipa_common$nr_genre %in% c("3")] <- "3 Genre"
ipa_common$nr_genreSum[ipa_common$nr_genre %in% c("4","5", "6","7","8","9")] <- "4 or more Genre"

# Genre-identification ----
# Generate frequency table 

tabyl(ipa_common, nr_genreSum, )
tabyl(ipa_common, nr_genreSum, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Mean, variance and standard deviation

describe(ipa_common$nr_genre)

# Comparison of countries and regions

tabyl(ipa_genre_all, genre_all, eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

tabyl(ipa_genre_all, genre_all, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

# Genre distribution ----
# Generate frequency table 

table(ipa_genre_all$genre_all)
prozent_genre <- round(table(ipa_genre_all$genre_all) / sum(table(ipa_genre_all$genre_all)) * 100, 2)
prozent_genre
sort(prozent_genre, decreasing = TRUE) 

# clean labels

ipa_genre_all <- ipa_genre_all %>% 
  mutate(
    genre_all = factor(genre_all,
                 levels = c("theater", 
                            "performance", 
                            "other_genre",
                            "music_theater", 
                            "interdis", 
                            "dance",
                            "comercial", 
                            "circus", 
                            "child_youth"),
                 labels = c("Theater",
                            "Performance",
                            "Other genre",
                            "Music theater", 
                            "Interdisciplinary artistic work",
                            "Dance",
                            "Commercial artistic work", 
                            "Circus", 
                            "Children's and youth theater")))
    
# Generate XTabs plot   
  
PlotXTabs2(ipa_genre_all, 
           eaipa_country, 
           genre_all, 
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",

           legend.title = "Genre",
           legend.position = "top",
           bar.outline.color = "darkgrey",

           title = "",
           caption = "The count refers to the number of mentions. 
           Individual respondents may have mentioned more than one genre.")

# Goodness of fit test ----
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_genre_all, eaipa_country == "Austria")
obs <- table(Austria$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826 )
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.3603 --> in etwa gleich

Bulgaria <- subset(ipa_genre_all, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826 )
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.3548 --> in etwa gleich

CzechRepublic <- subset(ipa_genre_all, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826 )
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.01243 --> 0Hyp verwerfen.. Weicht ab

Finland <- subset(ipa_genre_all, eaipa_country == "Finland")
obs <- table(Finland$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826 )
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.09521 --> in etwa gleich

Germany <- subset(ipa_genre_all, eaipa_country == "Germany")
obs <- table(Germany$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826 )
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 5.902e-05--> 0Hyp verwerfen.. Weicht ab

Hungary <- subset(ipa_genre_all, eaipa_country == "Hungary")
obs <- table(Hungary$genre_all)
exp <- c(0.1172,	0.0449,	0.1254,	0.1213,	0.0617,	0.0579,	0.1855,	0.2861)
chisq.test(x=obs, p=exp)
Hungary %>% count(genre_all)
# commercial fehlt - aber wenn excluded dann p=0.151

Iceland <- subset(ipa_genre_all, eaipa_country == "Iceland")
obs <- table(Iceland$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.003906 --> 0Hyp verwerfen.. Weicht ab

Italy <- subset(ipa_genre_all, eaipa_country == "Italy")
obs <- table(Italy$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.01821 --> 0Hyp verwerfen.. Weicht ab

Romania <- subset(ipa_genre_all, eaipa_country == "Romania")
obs <- table(Romania$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 2.328e-06 --> 0Hyp verwerfen.. Weicht ab

Slovenia <- subset(ipa_genre_all, eaipa_country == "Slovenia")
obs <- table(Slovenia$genre_all)
exp <- c(0.118875,	0.033175,	0.127075,	0.122975,	0.063375,	0.059575,	0.187175,	0.287775)
chisq.test(x=obs, p=exp)
Slovenia %>% count(genre_all)
# Circus is missing
# 0 hypothese in etwa gleich: p-value = 0.7057 --> 0Hyp verwerfen.. Weicht ab

Sweden <- subset(ipa_genre_all, eaipa_country == "Sweden")
obs <- table(Sweden$genre_all)
exp <- c(0.1172,	0.0449,	0.1254,	0.1213,	0.0617,	0.0579,	0.1855,	0.2861)
chisq.test(x=obs, p=exp)
Sweden %>% count(genre_all)
# commercial fehlt, aber wenn excluded dann: 0.151

Switzerland <- subset(ipa_genre_all, eaipa_country == "Switzerland")
obs <- table(Switzerland$genre_all)
exp <- c(0.1137, 0.0414, 0.0280, 0.1219, 0.1178, 0.0582, 0.0544, 0.1820, 0.2826)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.9041 --> in etwa gleich

# Profession ----

# Profession group allocation
ipa_ind$nr_prof <- as.integer(ipa_ind$nr_prof)
ipa_ind$nr_profSum[ipa_ind$nr_prof %in% c("1","0")] <- "1 profession"
ipa_ind$nr_profSum[ipa_ind$nr_prof %in% c("2")] <- "2 professions"
ipa_ind$nr_profSum[ipa_ind$nr_prof %in% c("3","4","5","6","7","8")] <- "3 or more professions"

# Generate frequency table 

tabyl(ipa_ind, nr_profSum, )
tabyl(ipa_ind, nr_profSum, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Mean, variance and standard deviation
describe(ipa_ind$nr_prof)

# Comparison of countries and regions

tabyl(ipa_profession, profession, eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

tabyl(ipa_profession, profession, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

# Goodness of fit test ----
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_profession, eaipa_country == "Austria")
obs <- table(Austria$profession)
exp <- c(0.0761,	0.4864,	0.0440,	0.0374,	0.0715,	0.1423,	0.1423)
chisq.test(x=obs, p=exp)
Austria %>% count(profession)
# techniciamn fehlt 
# 0 hypothese in etwa gleich: p-value = 0.4584 --> in etwa gleich

Bulgaria <- subset(ipa_profession, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$profession)
exp <- c(0.0728,	0.483,	0.0406,	0.034,	0.0681,	0.1389,	0.1389, 0.0237)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.1397 --> in etwa gleich

CzechRepublic <- subset(ipa_profession, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$profession)
exp <- c(0.0824,	0.4926,	0.0502,	0.0777,	0.1485,	0.1486)
chisq.test(x=obs, p=exp)
CzechRepublic %>% count(profession)
# designer und Technician fehlen. 
# 0 hypothese in etwa gleich: p-value = 0.0263 --> 0Hyp verwerfen.. Weicht ab

Finland <- subset(ipa_profession, eaipa_country == "Finland")
obs <- table(Finland$profession)
exp <- c(0.082528571,	0.492728571,	0.050328571,	0.043728571,	0.148628571,	0.148628571,	0.033428571)
chisq.test(x=obs, p=exp)
Finland %>% count(profession)
# Darmaturgie fehlt
# 0 hypothese in etwa gleich: p-value = 0.444 --> in etwa gleich

Germany <- subset(ipa_profession, eaipa_country == "Germany")
obs <- table(Germany$profession)
exp <- c(0.0728,	0.483,	0.0406,	0.034,	0.0681,	0.1389,	0.1389, 0.0237)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.2199 --> in etwa gleich

Hungary <- subset(ipa_profession, eaipa_country == "Hungary")
obs <- table(Hungary$profession)
exp <- c(0.1033,	0.5135,	0.1596,	0.1694,	0.0542)
chisq.test(x=obs, p=exp)
Hungary %>% count(profession)
# curator, designer, dramaturge fehlen
# commercial fehlt - aber wenn excluded dann p=0.1603 --> in etwa gleich

Iceland <- subset(ipa_profession, eaipa_country == "Iceland")
obs <- table(Iceland$profession)
exp <- c(0.0728,	0.483,	0.0406,	0.034,	0.0681,	0.1389,	0.1389, 0.0237)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.0002371 --> 0Hyp verwerfen.. Weicht ab

Italy <- subset(ipa_profession, eaipa_country == "Italy")
obs <- table(Italy$profession)
exp <- c(0.0728,	0.483,	0.0406,	0.034,	0.0681,	0.1389,	0.1389, 0.0237)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.02606 --> 0Hyp verwerfen.. Weicht ab

Romania <- subset(ipa_profession, eaipa_country == "Romania")
obs <- table(Romania$profession)
exp <- c(0.085233333,	0.495433333,	0.080533333,	0.151333333,	0.151333333,	0.036133333)
chisq.test(x=obs, p=exp)
Romania %>% count(profession)
# Curator & Designer fehlen
# 0 hypothese in etwa gleich: p-value = 0.1233 --> in etwa gleich

Slovenia <- subset(ipa_profession, eaipa_country == "Slovenia")
obs <- table(Slovenia$profession)
exp <- c(0.0728,	0.483,	0.0406,	0.034,	0.0681,	0.1389,	0.1389, 0.0237)
chisq.test(x=obs, p=exp)
Slovenia %>% count(profession)
# 0 hypothese in etwa gleich: p-value = 0.4742 --> in etwa gleich

Sweden <- subset(ipa_profession, eaipa_country == "Sweden")
obs <- table(Sweden$profession)
exp <- c(0.51248,	0.09758,	0.16838,	0.16838,	0.05318)
chisq.test(x=obs, p=exp)
Sweden %>% count(profession)
# designer curator artmanager fehlen
# commercial fehlt, aber wenn excluded dann: 0.1295

Switzerland <- subset(ipa_profession, eaipa_country == "Switzerland")
obs <- table(Switzerland$profession)
exp <- c(0.0728,	0.483,	0.0406,	0.034,	0.0681,	0.1389,	0.1389, 0.0237)
chisq.test(x=obs, p=exp)

#  Orgatype ----

# Generate frequency table 

tabyl(ipa_orgatype, nr_orgatype, )
tabyl(ipa_orgatype, nr_orgatype, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Mean, variance and standard deviation

describe(ipa_orgatype$nr_orgatype)

# Comparison of countries and regions

tabyl(ipa_orgatype, orgatype,eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

tabyl(ipa_orgatype, orgatype, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

# Generate XTabs plot   

PlotXTabs2(ipa_orgatype, 
           eaipa_country, 
           orgatype, 
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",
           
           legend.title = "Distribution of Types of Organization",
           legend.position = "top",
           bar.outline.color = "darkgrey",
           
           title = "",
           caption = "The count refers to the number of mentions. 
           Individual respondents may have mentioned more than one genre.")

# Goodness of fit test ----
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_orgatype, eaipa_country == "Austria")
obs <- table(Austria$orgatype)
exp <- c(0.49204,	0.12134,		0.10564,	0.07594,	0.08114,	0.06194,	0.06196)
chisq.test(x=obs, p=exp)
Austria %>% count(orgatype)
# orgatype other fehlt  
# 0 hypothese in etwa gleich: p-value = 0.9626 --> in etwa gleich

Bulgaria <- subset(ipa_orgatype, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.1468 --> in etwa gleich

CzechRepublic <- subset(ipa_orgatype, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
CzechRepublic %>% count(orgatype)
# 0 hypothese in etwa gleich: p-value = 0.002964 --> 0Hyp verwerfen.. Weicht ab

Finland <- subset(ipa_orgatype, eaipa_country == "Finland")
obs <- table(Finland$orgatype)
exp <- c(0.492042857,	0.121342857,	0.061942857,	0.105642857,	0.075942857,	0.081142857,		0.061942858)
chisq.test(x=obs, p=exp)
Finland %>% count(orgatype)
# venue_w_ensam fehlt
# 0 hypothese in etwa gleich: p-value = 0.9921 --> in etwa gleich

Germany <- subset(ipa_orgatype, eaipa_country == "Germany")
obs <- table(Germany$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.5501 --> in etwa gleich

Hungary <- subset(ipa_orgatype, eaipa_country == "Hungary")
obs <- table(Hungary$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
# p=0.6655 --> in etwa gleich

Iceland <- subset(ipa_orgatype, eaipa_country == "Iceland")
obs <- table(Iceland$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.9489 --> in etwa gleich

Italy <- subset(ipa_orgatype, eaipa_country == "Italy")
obs <- table(Italy$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.02693 --> 0Hyp verwerfen.. Weicht ab

Romania <- subset(ipa_orgatype, eaipa_country == "Romania")
obs <- table(Romania$orgatype)
exp <- c(0.4843,	0.1136,	0.0542,	0.0979,	0.0682,	0.0734,	0.0542,	0.0542)
chisq.test(x=obs, p=exp)
Romania %>% count(orgatype)
# 0 hypothese in etwa gleich: p-value = 0.9386 --> in etwa gleich

Slovenia <- subset(ipa_orgatype, eaipa_country == "Slovenia")
obs <- table(Slovenia$orgatype)
exp <- c(0.492042857,	0.121342857,	0.061942857,	0.105642857,	0.075942857,	0.081142857,		0.061942858)
chisq.test(x=obs, p=exp)
Slovenia %>% count(orgatype)
# venue_w_ensam fehlt
# 0 hypothese in etwa gleich: p-value = 0.2314 --> in etwa gleich

Sweden <- subset(ipa_orgatype, eaipa_country == "Sweden")
obs <- table(Sweden$orgatype)
exp <- c(0.500528571,	0.070428571,	0.114128571,	0.084428571,	0.089628571,	0.070428571,	0.070428574)
chisq.test(x=obs, p=exp)
Sweden %>% count(orgatype)
# dFestival fehlen
# commercial fehlt, aber wenn excluded dann: 0.1675

Switzerland <- subset(ipa_orgatype, eaipa_country == "Switzerland")
obs <- table(Switzerland$orgatype)
exp <- c(0.55335,	0.16695,	0.13725,	0.14245)
chisq.test(x=obs, p=exp)
Switzerland %>% count(orgatype)

# Legal status ---- 

# Filtering out artist companies

ipa_org_comp <- filter(ipa_org, company == 2)

# Generate frequency table 

prozent_legal_stat <- round(table(ipa_org_comp$legal_stat) / sum(table(ipa_org_comp$legal_stat)) * 100, 2)
prozent_legal_stat 
sort(prozent_legal_stat, decreasing = TRUE) 

# Comparison of countries and regions

tabyl(ipa_org_comp, legal_stat, eu_regions)
tabyl(ipa_org_comp, legal_stat, eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

tabyl(ipa_org_comp, legal_stat, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

# Generate XTabs plot   

PlotXTabs2(ipa_org_comp, 
           eaipa_country, 
           legal_stat, 
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",
           
           legend.title = "Distribution of Types of Legal Status",
           legend.position = "top",
           bar.outline.color = "darkgrey",
           
           title = "")

# Goodness of fit test ----
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_orgatype, eaipa_country == "Austria")
obs <- table(Austria$legal_stat)
exp <- c(0.5784667, 0.1680667, 0.2534666)
chisq.test(x=obs, p=exp)
Austria %>% count(legal_stat)
# Cooperative	Foundation	Limited liability company fehlt  
# 0 hypothese in etwa gleich: p-value = 0.007625 --> 0Hyp verwerfen.. Weicht ab

Bulgaria <- subset(ipa_orgatype, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$legal_stat)
exp <- c(0.5041,	0.0771,	0.0496,	0.0964,	0.0937,	0.1791)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.001976 --> 0Hyp verwerfen.. Weicht ab

CzechRepublic <- subset(ipa_orgatype, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$legal_stat)
exp <- c(0.5041,	0.0771,	0.0496,	0.0964,	0.0937,	0.1791)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 1.206e-12--> 0Hyp verwerfen.. Weicht ab

Finland <- subset(ipa_orgatype, eaipa_country == "Finland")
obs <- table(Finland$legal_stat)
exp <- c(0.56815, 0.11365, 0.16045, 0.15775)
chisq.test(x=obs, p=exp)
Finland %>% count(legal_stat)
# Cooperative and  fehlt
# 0 hypothese in etwa gleich: p-value = 0.04298--> 0Hyp verwerfen.. Weicht ab

Germany <- subset(ipa_orgatype, eaipa_country == "Germany")
obs <- table(Germany$legal_stat)
exp <- c(0.5041,	0.0771,	0.0496,	0.0964,	0.0937,	0.1791)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.001247 --> 0Hyp verwerfen.. Weicht ab

Hungary <- subset(ipa_orgatype, eaipa_country == "Hungary")
obs <- table(Hungary$legal_stat)
exp <- c(0.5468, 0.0923, 0.1391, 0.2218)
chisq.test(x=obs, p=exp)
Hungary %>% count(legal_stat)
# other und cooperative fehlen
# p=1.738e-05 --> 0Hyp verwerfen.. Weicht ab

Iceland <- subset(ipa_orgatype, eaipa_country == "Iceland")
obs <- table(Iceland$legal_stat)
exp <- c(0.5041,	0.0771,	0.0496,	0.0964,	0.0937,	0.1791)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.1351 --> in etwa gleich

Italy <- subset(ipa_orgatype, eaipa_country == "Italy")
obs <- table(Italy$legal_stat)
exp <- c(0.5041,	0.0771,	0.0496,	0.0964,	0.0937,	0.1791)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 1.964e-09 --> 0Hyp verwerfen.. Weicht ab

Romania <- subset(ipa_orgatype, eaipa_country == "Romania")
obs <- table(Romania$legal_stat)
exp <- c(0.6625, 0.3375)
chisq.test(x=obs, p=exp)
Romania %>% count(legal_stat)
# Cooperative	Foundation	Limited liability company	Other fehlt
# 0 hypothese in etwa gleich: p-value = 0.301 --> in etwa gleich

Slovenia <- subset(ipa_orgatype, eaipa_country == "Slovenia")
obs <- table(Slovenia$legal_stat)
exp <- c(0.6625, 0.3375)
chisq.test(x=obs, p=exp)
Slovenia %>% count(legal_stat)
# Cooperative	Foundation	Limited liability company	Other fehlt
# 0 hypothese in etwa gleich: p-value = 0.002549 --> 0Hyp verwerfen.. Weicht ab

Sweden <- subset(ipa_orgatype, eaipa_country == "Sweden")
obs <- table(Sweden$legal_stat)
exp <- c(0.5041,	0.0771,	0.0496,	0.0964,	0.0937,	0.1791)
chisq.test(x=obs, p=exp)
Sweden %>% count(legal_stat)
# commercial fehlt, aber wenn excluded dann: 6.239e-05 --> 0Hyp verwerfen.. Weicht ab

Switzerland <- subset(ipa_orgatype, eaipa_country == "Switzerland")
obs <- table(Switzerland$legal_stat)
exp <- c(0.535775, 0.128075,	0.125375,	0.210775)
chisq.test(x=obs, p=exp)
Switzerland %>% count(legal_stat)
# Cooperative und foundation fehlen
# 0 hypothese in etwa gleich: p-value = 0.6277 --> in etwa gleich

# NPO status -----

# Generate frequency table 

table(ipa_org$npo_stat)
prozent<- round(table(ipa_org$npo_stat) / sum(table(ipa_org$npo_stat)) * 100, 2)
prozent
sort(prozent, decreasing = TRUE) 

# Comparison of countries and regions

table(ipa_org$eaipa_country, ipa_org$npo_stat)
tabyl(ipa_org, npo_stat, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Generate XTabs plot   

PlotXTabs2(ipa_orgatype, 
           eaipa_country, 
           npo_stat, 
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",
           
           legend.title = "NPO Status",
           legend.position = "top",
           bar.outline.color = "darkgrey",
           
           title = "",)

# Goodness of fit test ----
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_org, eaipa_country == "Austria")
obs <- table(Austria$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.05672 --> Gleich

Bulgaria <- subset(ipa_org, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.5649 --> Gleich

CzechRepublic <- subset(ipa_org, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.6565--> Gleich 

Finland <- subset(ipa_org, eaipa_country == "Finland")
obs <- table(Finland$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.7227--> Gleich

Germany <- subset(ipa_org, eaipa_country == "Germany")
obs <- table(Germany$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 2.521e-08 --> 0Hyp verwerfen.. Weicht ab

Hungary <- subset(ipa_org, eaipa_country == "Hungary")
obs <- table(Hungary$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# p=0.06893 --> gleich 

Iceland <- subset(ipa_org, eaipa_country == "Iceland")
obs <- table(Iceland$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.4725 --> Gleich

Italy <- subset(ipa_org, eaipa_country == "Italy")
obs <- table(Italy$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.627 --> Gleich

Romania <- subset(ipa_org, eaipa_country == "Romania")
obs <- table(Romania$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.5221 --> Gleich

Slovenia <- subset(ipa_org, eaipa_country == "Slovenia")
obs <- table(Slovenia$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.002549 --> 0Hyp verwerfen.. Weicht ab

Sweden <- subset(ipa_org, eaipa_country == "Sweden")
obs <- table(Sweden$npo_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
# commercial fehlt, aber wenn excluded dann: 0.0491 --> 0Hyp verwerfen.. Weicht ab

Switzerland <- subset(ipa_org, eaipa_country == "Switzerland")
obs <- table(Switzerland$legal_stat)
exp <- c(0.1961, 0.8039)
chisq.test(x=obs, p=exp)
Switzerland %>% count(legal_stat)
# Cooperative und foundation fehlen
# 0 hypothese in etwa gleich: p-value = 0.6277 --> in etwa gleich

# Employment type ----

# Generate frequency table 

prozent_employment <- round(table(ipa_ind$employment) / sum(table(ipa_ind$employment)) * 100, 2)
prozent_employment
sort(prozent_employment, decreasing = TRUE) 

# Comparison of countries and regions

tabyl(ipa_ind, employment, eu_regions)
tabyl(ipa_ind, employment, eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

tabyl(ipa_ind, employment, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) %>%
  adorn_pct_formatting(digits = 1) 

# Goodness of fit test 
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_ind, eaipa_country == "Austria")
obs <- table(Austria$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.009929 --> weicht ab

Bulgaria <- subset(ipa_ind, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.6717 --> Gleich

CzechRepublic <- subset(ipa_ind, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$employment)
exp <- c(0.652133333,	0.100233333,	0.247633334)
chisq.test(x=obs, p=exp)
CzechRepublic %>% count(employment)
# part time fehlt
# 0 hypothese in etwa gleich: p-value = 0.6598--> Gleich 

Finland <- subset(ipa_ind, eaipa_country == "Finland")
obs <- table(Finland$employment)
exp <- c(0.664166667,	0.259666666, 0.076166667)
chisq.test(x=obs, p=exp)
Finland %>% count(employment)
# full time fehlt
# 0 hypothese in etwa gleich: p-value = 0.6274--> Gleich

Germany <- subset(ipa_ind, eaipa_country == "Germany")
obs <- table(Germany$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
Germany %>% count(employment)
# 0 hypothese in etwa gleich: p-value = 0.001996 --> 0Hyp verwerfen.. Weicht ab

Hungary <- subset(ipa_ind, eaipa_country == "Hungary")
obs <- table(Hungary$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# p=0.06473 --> gleich 

Iceland <- subset(ipa_ind, eaipa_country == "Iceland")
obs <- table(Iceland$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.9619 --> Gleich

Italy <- subset(ipa_ind, eaipa_country == "Italy")
obs <- table(Italy$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 1.421e-06 Weicht ab

Romania <- subset(ipa_ind, eaipa_country == "Romania")
obs <- table(Romania$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.808 --> Gleich

Slovenia <- subset(ipa_ind, eaipa_country == "Slovenia")
obs <- table(Slovenia$employment)
exp <- c(0.652133333,	0.100233333,	0.247633334)
chisq.test(x=obs, p=exp)
Slovenia %>% count(employment)
# parttime fehlt
# 0 hypothese in etwa gleich: p-value = 0.2657 --> gleich

Sweden <- subset(ipa_ind, eaipa_country == "Sweden")
obs <- table(Sweden$employment)
exp <- c(0.652133333,	0.100233333,	0.247633334)
chisq.test(x=obs, p=exp)
Sweden %>% count(employment)
# Partime fehlt
# commercial fehlt, aber wenn excluded dann: 0.026 --> 0Hyp verwerfen.. Weicht ab

Switzerland <- subset(ipa_ind, eaipa_country == "Switzerland")
obs <- table(Switzerland$employment)
exp <- c(0.6361, 0.0842, 0.2316, 0.0481)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.05647 --> gleich

#Employment ----

# Share of employment forms in organizations

prozent <- round(table(ipa_org$full_time19) / sum(table(ipa_org$full_time19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$part_time19) / sum(table(ipa_org$part_time19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$temporary19) / sum(table(ipa_org$temporary19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$fee19) / sum(table(ipa_org$fee19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$special_order19) / sum(table(ipa_org$special_order19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$volun_intern19) / sum(table(ipa_org$special_order19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$volun_intern19) / sum(table(ipa_org$special_order19)) * 100, 2)
prozent
prozent<- round(table(ipa_org$no_cot19) / sum(table(ipa_org$no_cot19)) * 100, 2)
prozent

# Location ------

# Generate frequency table 

table(ipa_common$location)
prozent_location <- round(table(ipa_common$location) / sum(table(ipa_common$location)) * 100, 2)
prozent_location
sort(prozent_location, decreasing = TRUE) 

# Comparison of countries and regions

table(ipa_common$location, ipa_common$eaipa_country)
tabyl(ipa_common, location, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1) 

# Pie chart location distribution

prozent_nr <- round(table(ipa_common$location) / sum(table(ipa_common$location)) * 100, 2)
Beschriftung_nr <- c("City", "Large city", "Rural area", "Town")
Beschriftung_nr <- paste(Beschriftung_nr, prozent_nr, "%", " ")

labels(ipa_common$location)
pie(table(ipa_common$location), labels = Beschriftung_nr, col=gray(seq(0.3, 1.0, length = 4)))

# clean labels

ipa_common <- ipa_common %>% 
  mutate(
    location = factor(location,
                       levels = c("town", 
                                  "rural", 
                                  "large_city",
                                  "city"
                                  ),
                       labels = c("Town",
                                  "Rural area",
                                  "Large city",
                                  "City")))

# Generate XTabs plot

PlotXTabs2(ipa_common, 
           eaipa_country, 
           location, 
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",
           
           legend.title = "Location",
           legend.position = "top",
           bar.outline.color = "darkgrey",
           
           title = "",)

# Goodness of fit test 
# Adjustment of the distributions in case of missing level

Austria <- subset(ipa_common, eaipa_country == "Austria")
obs <- table(Austria$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.02298 --> weicht ab

Bulgaria <- subset(ipa_common, eaipa_country == "Bulgaria")
obs <- table(Bulgaria$location)
exp <- c(0.145533333,	0.788833333,		0.065633333)
chisq.test(x=obs, p=exp)
Bulgaria %>% count(location)
# Rural fehlt
# 0 hypothese in etwa gleich: p-value = 0.02826 --> weicht ab

CzechRepublic <- subset(ipa_common, eaipa_country == "Czech Republic")
obs <- table(CzechRepublic$location)
exp <- c(0.140966667,	0.784266667,	0.074766667)
chisq.test(x=obs, p=exp)
CzechRepublic %>% count(location)
# Town fehlt
# 0 hypothese in etwa gleich: p-value = 0.03357--> weicht ab 

Finland <- subset(ipa_common, eaipa_country == "Finland")
obs <- table(Finland$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
Finland %>% count(location)
# 0 hypothese in etwa gleich: p-value = 0.489--> Gleich

Germany <- subset(ipa_common, eaipa_country == "Germany")
obs <- table(Germany$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
Germany %>% count(location)
# 0 hypothese in etwa gleich: p-value = 0.1757 --> gleich

Hungary <- subset(ipa_common, eaipa_country == "Hungary")
obs <- table(Hungary$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
# p=0.1757 --> gleich 

Iceland <- subset(ipa_common, eaipa_country == "Iceland")
obs <- table(Iceland$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 0.8341 --> Gleich

Italy <- subset(ipa_common, eaipa_country == "Italy")
obs <- table(Italy$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 3.871e-09 Weicht ab

Romania <- subset(ipa_common, eaipa_country == "Romania")
obs <- table(Romania$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
Romania %>% count(location)
# Rural und Town fehlen

Slovenia <- subset(ipa_common, eaipa_country == "Slovenia")
obs <- table(Slovenia$location)
exp <- c(0.140966667,	0.784266667,	0.074766667)
chisq.test(x=obs, p=exp)
Slovenia %>% count(location)
# Town fehlt
# 0 hypothese in etwa gleich: p-value = 0.2737 --> gleich

Sweden <- subset(ipa_common, eaipa_country == "Sweden")
obs <- table(Sweden$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
Sweden %>% count(location)
# commercial fehlt, aber wenn excluded dann: 0.4964 --> Gleich

Switzerland <- subset(ipa_common, eaipa_country == "Switzerland")
obs <- table(Switzerland$location)
exp <- c(0.1257,	0.769,	0.0595,	0.0458)
chisq.test(x=obs, p=exp)
# 0 hypothese in etwa gleich: p-value = 1.702e-05--> weicht ab

# Archetype Performing Arts production ----
# Modus von Location, emplyment, & Legal status

Mode(na.omit(ipa_common$location))
Mode(na.omit(ipa_ind$employment))
Mode(na.omit(ipa_org$legal_stat))
Mode(na.omit(ipa_orgatype$orgatype))
Mode(na.omit(ipa_profession$profession))

# Instuitutional definition ----

# Describe individual items

describe(ipa_common$incr_info)
describe(ipa_common$need_f_com)
describe(ipa_common$need_f_collab)
describe(ipa_common$pol_alliance)
describe(ipa_common$build_struc)
describe(ipa_common$joint_venture)

# Build structuration scale 

ipa_struc <- subset(ipa_common, select = c(case, need_f_com:joint_venture, particip_i_advo))

str(ipa_struc)
summary(ipa_struc)

# % NA per row

percentmissing = function(x){ sum(is.na(x))/length(x)*100}
missing = apply(ipa_struc, 1, percentmissing)
table(missing)

# delete rows with more than 20% NA

replaceorg = subset(ipa_struc, missing <= 20)

# % NA in column
apply(replaceorg, 2, percentmissing)

# Multiple impute datapoints 
tempnomiss = mice(replaceorg) 
nomiss = complete(tempnomiss, 1)
summary(nomiss)

structuration = c(2, 3, 4, 5, 6, 7, 8, 9) 

# Estimates of reliability

psych::alpha(ipa_struc[structuration])

#create new factor score

ipa_struc$structuration = apply(ipa_struc[ ,structuration], 1, mean)

# merge Datasets according to case 

ipa_common_final <- merge( ipa_common, ipa_struc, by = "case")

# Describe structuration index

describe(ipa_common_final$structuration)

# Compare countries and regions 

Beschriftung_nr <- c("Identification with one genre", " Identification with two genre", "Identification with three Genre", "Identification with four or more Genre")
Beschriftung_nr <- paste(Beschriftung_nr, prozent_nr, "%", " ")
pie(table(ipa_common$nr_genreSum), labels = Beschriftung_nr, main = "Identification with multiple Genres")

boxplot(ipa_common_final$structuration ~ ipa_common_final$eu_regions, 
        xlab = "", ylab = "Structuration level")

describeBy(ipa_common_final$structuration, ipa_common_final$eu_regions)

# Describe structuration for regions

struEU <- subset(ipa_common_final, eu_regions == "Eastern Europe")
median(struEU$structuration, na.rm = TRUE)
struEU <- subset(ipa_common_final, eu_regions == "Northern Europe")
median(struEU$structuration, na.rm = TRUE)
struEU <- subset(ipa_common_final, eu_regions == "Southern Europe")
median(struEU$structuration, na.rm = TRUE)
struEU <- subset(ipa_common_final, eu_regions == "Western Europe")
median(struEU$structuration, na.rm = TRUE)

# Output ------
# Output per artistic company and per facilitating organization in country comparison

# Premieres
ipa_org_comp <- filter(ipa_org, company == 2)
describe(ipa_org_comp$prem19)
describeBy(ipa_org_comp$prem19, group = ipa_org_comp$eaipa_country)
ipa_org_org <- filter(ipa_org, product_house == 2 | rehear_space == 2 | venue_w_ensam == 2 | venue_wo_ensam == 2 | festival == 2 | resid_space == 2)
describe(ipa_org_org$prem19)
describeBy(ipa_org_org$prem19, group = ipa_org_org$eu_regions, skew = FALSE)

# Shows
ipa_org_comp <- filter(ipa_org, company == 2)
describe(ipa_org_comp$shows19)
describeBy(ipa_org_comp$shows19, group = ipa_org_comp$eu_regions)
ipa_org_org <- filter(ipa_org, product_house == 2 | rehear_space == 2 | venue_w_ensam == 2 | venue_wo_ensam == 2 | festival == 2 | resid_space == 2)
describe(ipa_org_org$shows19)
describeBy(ipa_org_org$shows19, group = ipa_org_org$eaipa_country)

# spectators
ipa_org_comp <- filter(ipa_org, company == 2)
describe(ipa_org_comp$spectators19)
describeBy(ipa_org_comp$spectators19, group = ipa_org_comp$eu_regions)
ipa_org_org <- filter(ipa_org, product_house == 2 | rehear_space == 2 | venue_w_ensam == 2 | venue_wo_ensam == 2 | festival == 2 | resid_space == 2)
describe(ipa_org_org$spectators19)
describeBy(ipa_org_org$spectators19, group = ipa_org_org$eu_regions)

# Shows_avre19
ipa_org_comp <- filter(ipa_org, company == 2)
describe(ipa_org_comp$shows_avre19)
describeBy(ipa_org_comp$shows_avre19, group = ipa_org_comp$eu_regions)
ipa_org_org <- filter(ipa_org, product_house == 2 | rehear_space == 2 | venue_w_ensam == 2 | venue_wo_ensam == 2 | festival == 2 | resid_space == 2)
describe(ipa_org_org$shows_avre19)
describeBy(ipa_org_org$shows_avre19, group = ipa_org_org$eu_regions, skew = FALSE)

# Crrelation Structuration / output-------
# um es zu behandeln als wären es Intervall daten muss 0 = 0 sein 

ipa_org_comp$output <- rowMeans(subset(ipa_org_comp, select = c(spectators19, shows19,prem19, shows_avre19)), na.rm = TRUE)

# Correlation structuration and output

cor.test(ipa_org$structuration, ipa_org$shows_avre19, method = "spearman")
cor.test(ipa_org$structuration, ipa_org$guestperf19, method = "spearman")

# Funding ----

# Read in data 
exp_data <- read_excel("~/Documents/_PhD Promotion/Analyses/Qual Analyse/Daten/data_Interview_IPA_Association_2021-04-19_13-25.xlsx")
ipa_common <- read_csv("Data/Data_out/ipa_common.csv")
ipa_org <- read_csv("Data/Data_out/ipa_org.csv")
ipa_ind <- read_csv("Data/Data_out/ipa_ind.csv")
ipa_orgatype <- read_csv("Data/Data_out/_prof_org/ipa_orgatype.csv")

# Clean data from the expert survey
 
exp_data <- rename(exp_data,
                   proj_fun_nat = RA06_01,
                   schol_nat = RA06_02,
                   struct_fun_nat = RA06_03,
                   newc_fun_nat = RA06_04, 
                   mobi_nat = RA06_05,
                   int_mobi_nat = RA06_06,
                   venue_fund_nat = RA06_07,
                   proj_fun_mun = RA07_01,
                   schol_mun = RA07_02,
                   struct_fun_mun = RA07_03,
                   newc_fun_mun = RA07_04, 
                   mobi_mun = RA07_05,
                   int_mobi_mun = RA07_06,
                   venue_fund_mun = RA07_07)

exp_data <- clean_names(exp_data)
exp_data <- exp_data[-c(8,9,13), ]

exp_data$eu_regions[exp_data$a002_02 %in% c("Austria",
                                              "Germany",
                                              "Switzerland")] <- "Western Europe"

exp_data$eu_regions[exp_data$a002_02 %in% c("Bulgaria",
                                              "Czech Republic",
                                              "Hungary",
                                              "Romania"
                                              )] <- "Former Sovjet Countries"

exp_data$eu_regions[exp_data$a002_02 %in% c("Italy",
                                            "Slovenia"
                                            )] <- "Southwestern Europe"

exp_data$eu_regions[exp_data$a002_02 %in% c("Finland",
                                              "Iceland",
                                              "Sweden")] <- "Scandinavia"

exp_data_fun_nat <- exp_data %>% pivot_longer(names_to = "fun_nat", values_to = "yes_no1", proj_fun_nat:venue_fund_nat)
exp_data_fun_nat <- filter(exp_data_fun_nat, yes_no1 == 1)

exp_data_fun_mun <- exp_data %>% pivot_longer(names_to = "fun_mun", values_to = "yes_no2", proj_fun_mun:venue_fund_mun)
exp_data_fun_mun <- filter(exp_data_fun_mun, yes_no2 == 1)

# Differences in funding instruments

table(exp_data_fun_nat$a002_02, exp_data_fun_nat$fun_nat)
table(exp_data_fun_mun$a002_02, exp_data_fun_mun$fun_mun)

# Frequencies of the funding instruments

barplot(sort(table(exp_data_fun_nat$fun_nat), decreasing = TRUE), names.arg=c("Project", "Mobility int.", "New talent", "Mobility nat.", "Venue", "Scholarships", "Structure"), xlab = "", ylab = "Number of countries")
barplot(sort(table(exp_data_fun_mun$fun_mun), decreasing = TRUE), names.arg=c("Project", "Venue", "Mobility int.", "Structure", "New talent", "Scholarships", "Mobility nat."), xlab = "", ylab = "Number of capitals")

# Public funding in organization budget

# clean labels 
ipa_org$fund[ipa_org$pub_fund_org %in% c("1")] <- "90 - 100 %"
ipa_org$fund[ipa_org$pub_fund_org %in% c("2")] <- "50 - 89 %"
ipa_org$fund[ipa_org$pub_fund_org %in% c("3")] <- "20 - 49 %"
ipa_org$fund[ipa_org$pub_fund_org %in% c("4")] <- "1 - 20 %"
ipa_org$fund[ipa_org$pub_fund_org %in% c("5")] <- "0 %"

prozent<- round(table(ipa_org$fund) / sum(table(ipa_org$fund)) * 100, 2)
prozent

describe(ipa_org$pub_fund_org)

# country and regional differences

describeBy(ipa_org$pub_fund_org, group = ipa_org$eaipa_country)
describeBy(ipa_org$pub_fund_org, group = ipa_org$eu_regions)

table(ipa_org$eaipa_country, ipa_org$fund)
tabyl(ipa_org, fund, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Bar chart public funding in organization budget country

data <- ipa_org[, -2:-120] 
data <- na.omit(data)

ggplot(data = data, aes(y = eaipa_country, fill = fund)) +
        geom_bar()+
  labs(x="Number of Organizations", y = "Country", fill = "Share of public funding relative to total income")

# Bar chart public funding in organization budget EU region

data <- ipa_org[, -2:-120] 
data <- na.omit(data)

ggplot(data = data, aes(y = eu_regions, fill = fund)) +
  geom_bar()+
  labs(x="Number of Organizations", y = "EU Region", fill = "Share of public funding relative to total income")

#SOCIAL SITUATION

# Income ----

# Clean data

ipa_ind_money <- ipa_ind[ , c(1, 65:66)]

# % NA in rows

percentmissing = function(x){ sum(is.na(x))/length(x)*100}
missing = apply(ipa_ind_money, 1, percentmissing)
table(missing)

# Delete with more 35% NA

replacemoney = subset(ipa_ind_money, missing <= 35)

# % NA in coloms

apply(replacemoney, 2, percentmissing)

# Multiple impute datapoints 

tempnomiss = mice(replacemoney) 
nomiss = complete(tempnomiss, 1)
summary(nomiss)

#exclude outliers

cutoff = qchisq(1-.001, ncol(nomiss))
mahal = mahalanobis(nomiss,
                    colMeans(nomiss),
                    cov(nomiss))

cutoff # cutoff score

ncol(nomiss) # df
summary(mahal < cutoff)

# cut off outliers 

noout = subset(nomiss, mahal < cutoff)

# merge Dataset according to case

ipa_money <- merge(ipa_ind, noout, by = "case")

# Plot income according to country 

ggplot(ipa_money, aes(x=income19_eur.x, y=eaipa_country)) + 
  geom_boxplot()

# Plot income according to gender

# Clean labels 

ipa_money <- ipa_money %>% 
  mutate(
    sex = factor(sex,
                         levels = c("Female", "Male", "other"),
                         labels = c("Female","Male","Diverse")))

ggplot(data = ipa_money, mapping = aes(x=income19_eur.x, 
                      y=eaipa_country, fill = sex)) + 
  geom_boxplot() +
  theme_linedraw() + 
  theme(legend.position = "bottom") + 
  scale_fill_brewer(palette="Greys") +
  labs(title = "",
       x = "Individual net income per year in €",
       y = "",
       fill = "Gender",
       caption = "") 
  
# Plot income according to employment

# Clean labels

ipa_money <- ipa_money %>% 
  mutate(
    employment = factor(employment,
                 levels = c("freelancer", 
                            "full-time employee", 
                            "hybrid (freelance/empl)", 
                            "part-time employee"),
                 labels = c("Self-employed",
                            "Full time employed",
                            "Hybrid employed",
                            "Part time employed")))

ggplot(ipa_money, aes(x=income19_eur.x, y=eaipa_country, fill=employment)) + 
  geom_boxplot()

# Descriptive statistics
describe(ipa_money$income20_eur.x)
describeBy(ipa_money$income19_eur.x, ipa_money$eaipa_country)

# Share of professionals below the at-risk-of-poverty threshold ----

describe(ipa_money$income19_eur.x)
describe(ipa_money$income20_eur.x)

# Germany ----

ipa_germany <- filter(ipa_money, country == "Germany") %>%
  subset(select = c(income19_eur.x, income20_eur.x, sex))  %>%
  mutate(poor19 = income19_eur.x < 14109,
         poor20 = income20_eur.x < 14109)

# median income of IPA professionals

median(ipa_germany$income19_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_germany,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_germany,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_germany,poor19, sex ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_germany,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind, degree, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Austria ----

ipa_ind_austria <- filter(ipa_money, country == "Austria") %>%
  mutate(poor19 = income19_eur.x < 15437,
         poor20 = income20_eur.x < 15437)

# median income of IPA professionals

median(ipa_ind_austria$income19_eur.x, na.rm = TRUE)
median(ipa_ind_austria$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_austria,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_austria,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_austria,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_austria,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Iceland ----

ipa_ind_iceland <- filter(ipa_money, country == "Iceland") %>%
  mutate(poor19 = income19_eur.x < 23951,
         poor20 = income20_eur.x < 23951)

# median income of IPA professionals

median(ipa_ind_iceland$income19_eur.x, na.rm = TRUE)
median(ipa_ind_iceland$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_iceland,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_iceland,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_iceland,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_iceland,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Italy ----

ipa_ind_italy <- filter(ipa_money, country == "Italy") %>%
  mutate(poor19 = income19_eur.x < 10299,
         poor20 = income20_eur.x < 10299)

# median income of IPA professionals

median(ipa_ind_italy$income19_eur.x, na.rm = TRUE)
median(ipa_ind_italy$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_italy,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_italy,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_italy,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_italy,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Bulgaria ----

ipa_ind_bulgaria <- filter(ipa_money, country == "Bulgaria") %>%
  mutate(poor19 = income19_eur.x < 2534,
         poor20 = income20_eur.x < 2762)

# median income of IPA professionals

median(ipa_ind_bulgaria$income19_eur.x, na.rm = TRUE)
median(ipa_ind_bulgaria$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_bulgaria,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_bulgaria,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_bulgaria,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_bulgaria,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Switzerland ----

ipa_ind_switzerland <- filter(ipa_money, country == "Switzerland") %>%
  mutate(poor19 = income19_eur.x < 24755,
         poor20 = income20_eur.x < 24755)

# median income of IPA professionals

median(ipa_ind_switzerland$income19_eur.x, na.rm = TRUE)
median(ipa_ind_switzerland$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_switzerland,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_switzerland,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_switzerland,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_switzerland,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Romania ----

ipa_ind_romania <- filter(ipa_money, country == "Romania") %>%
  mutate(poor19 = income19_eur.x < 2310,
         poor20 = income20_eur.x < 2310)

# median income of IPA professionals

median(ipa_ind_romania$income19_eur.x, na.rm = TRUE)
median(ipa_ind_romania$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_romania,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_romania,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_romania,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_romania,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Hungary ----

ipa_ind_hungary <- filter(ipa_money, country == "Hungary") %>%
  mutate(poor19 = income19_eur.x < 3511,
         poor20 = income20_eur.x < 3511)

# median income of IPA professionals

median(ipa_ind_hungary$income19_eur.x, na.rm = TRUE)
median(ipa_ind_hungary$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_hungary,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_hungary,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_hungary,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_hungary,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Sweden ----

ipa_ind_sweden <- filter(ipa_money, country == "Sweden") %>%
  mutate(poor19 = income19_eur.x < 14684,
         poor20 = income20_eur.x < 14684)

# median income of IPA professionals

median(ipa_ind_sweden$income19_eur.x, na.rm = TRUE)
median(ipa_ind_sweden$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_sweden,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_sweden,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_sweden,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_sweden,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Finland ----

ipa_ind_finland <- filter(ipa_money, country == "Finland") %>%
  mutate(poor19 = income19_eur.x < 14927,
         poor20 = income20_eur.x < 14927)

# median income of IPA professionals

median(ipa_ind_finland$income19_eur.x, na.rm = TRUE)
median(ipa_ind_finland$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_finland,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_finland,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_finland,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_finland,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Czech Republic ----

ipa_ind_czech <- filter(ipa_money, country == "Czech Republic") %>%
  mutate(poor19 = income19_eur.x < 5997,
         poor20 = income20_eur.x < 5997)

# median income of IPA professionals

median(ipa_ind_czech$income19_eur.x, na.rm = TRUE)
median(ipa_ind_czech$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_czech,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_czech,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_czech,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_czech,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Slovenia ----

ipa_ind_slovenia <- filter(ipa_money, country == "Slovenia") %>%
  mutate(poor19 = income19_eur.x < 8440,
         poor20 = income20_eur.x < 8440)

# median income of IPA professionals

median(ipa_ind_slovenia$income19_eur.x, na.rm = TRUE)
median(ipa_ind_slovenia$income20_eur.x, na.rm = TRUE)

# % below the at-risk-of-poverty threshold

tabyl(ipa_ind_slovenia,poor19, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_slovenia,poor20, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_slovenia,poor19, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind_slovenia,poor20, sex) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Income from a second profession to secure  livelihood ----

ipa_ind <- read_csv("Data/Data_out/ipa_ind.csv")

# name values
ipa_ind$empl_ouside_ipa[ipa_ind$empl_ouside_ipa %in% c("1")] <- "0 %"
ipa_ind$empl_ouside_ipa[ipa_ind$empl_ouside_ipa %in% c("2")] <- "1 - 25 %"
ipa_ind$empl_ouside_ipa[ipa_ind$empl_ouside_ipa %in% c("3")] <- "26 - 50 %"
ipa_ind$empl_ouside_ipa[ipa_ind$empl_ouside_ipa %in% c("4")] <- "51 - 75 %"
ipa_ind$empl_ouside_ipa[ipa_ind$empl_ouside_ipa %in% c("5")] <- "75 - 99 %"
ipa_ind$empl_ouside_ipa[ipa_ind$empl_ouside_ipa %in% c("6")] <- "all"

# % of values

table(ipa_ind$empl_ouside_ipa)
prozent<- round(table(ipa_ind$empl_ouside_ipa) / sum(table(ipa_ind$empl_ouside_ipa)) * 100, 2)
prozent

# % of values per country

table(ipa_ind$eaipa_country, ipa_ind$empl_ouside_ipa)
tabyl(ipa_ind, empl_ouside_ipa, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Coverage by social security ----

# Not covered by social security (1= not selected, 2=selected)

# Compare regions 

tabyl(ipa_ind, no_soci_sec, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Compare countries 

tabyl(ipa_ind, no_soci_sec, eaipa_country)
tabyl(ipa_ind, no_soci_sec, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Clean names

ipa_ind <- ipa_ind %>% 
  mutate(
    no_soci_sec = factor(no_soci_sec,
                       levels = c("1", "2"),
                       labels = c("Fully covered by social security",
                                  "Not fully covered by social security")))

# Plot country comparison 

PlotXTabs2(ipa_ind, 
           eaipa_country, 
           no_soci_sec,
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",
           
           legend.title = "Self-report concearning social security coverage",
           legend.position = "bottom",
           bar.outline.color = "darkgrey",
           
           title = "",
           caption = "")

# Frequency table social security for artists 

tabyl(ipa_ind, artist_soci_sec, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Frequency table unemployment insurance

describe(ipa_ind$unempl_insur)
tabyl(ipa_ind, unempl_insur, eaipa_country)
tabyl(ipa_ind, unempl_insur, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# Clean labels 

ipa_ind <- ipa_ind %>% 
  mutate(
    unempl_insur = factor(unempl_insur,
                         levels = c("1", "2"),
                         labels = c("Covered by unemployment insurance",
                                    "Not covered by unemployment insurance")))

# Stabs plot unemployment insurance

PlotXTabs2(ipa_ind, 
           eaipa_country, 
           unempl_insur,
           
           xlab = "",
           ylab = NULL,
           x.axis.orientation = "slant",
           
           sample.size.label = TRUE,
           label.fill.alpha = .6,
           label.text.size = 3,
           
           plottype = "percent",
           results.subtitle = FALSE,
           palette = "Greys",
           
           legend.title = "Self-report concearning unemployment insurance coverage",
           legend.position = "bottom",
           bar.outline.color = "darkgrey",
           
           title = "",
           caption = "")

# Quality of work ----  

# Build quality of work scale 

ipa_ind <- read_csv("Data/Data_out/ipa_ind.csv")
ipa_work_qual <- subset(ipa_ind, select = c(case, 
                                            income_plenty,    
                                            income_adequacy,  
                                            fin_secu, 
                                            job_stabi,
                                            pension_secu,
                                            unpaid_work,
                                            pos_toward_overwork,
                                            ident_w_prof,
                                            meaningfullness_o_work))

#revers coloms

cols=c("unpaid_work", "pos_toward_overwork")
ipa_work_qual [ ,cols] = 6 - ipa_work_qual[ ,cols]

str(ipa_work_qual)
summary(ipa_work_qual)

# % NA in rows

percentmissing = function(x){ sum(is.na(x))/length(x)*100}
missing = apply(ipa_work_qual, 1, percentmissing)
table(missing)

# Delete rows with more then 20% NA

replaceorg = subset(ipa_work_qual, missing <= 20)

# &NA in column

apply(replaceorg, 2, percentmissing)

# Multiple impute datapoints 

tempnomiss = mice(replaceorg) 
nomiss = complete(tempnomiss, 1)
summary(nomiss)

work_qual = c(2, 3, 4, 5, 6, 7,8,9,10) 

# Reliability test 

psych::alpha(ipa_work_qual[work_qual])

#create new factor score

ipa_work_qual$work_qual = apply(ipa_work_qual[ ,work_qual], 1, mean)

# Merge Datasets 

ipa_ind <- merge(ipa_ind, ipa_work_qual, by = "case")

# Describe work qual index

describe(ipa_ind$work_qual)
describeBy(ipa_ind$work_qual, ipa_ind$eu_regions)
describeBy(ipa_ind$work_qual, ipa_ind$eaipa_country)

# Correlations with age and education level

cor.test(ipa_ind$work_qual, ipa_ind$age)
cor.test(ipa_ind$work_qual, ipa_ind$degree)

# Income indicator

ipa_ind$income <- rowMeans(subset(ipa_ind, select = c(income_plenty.y,    
                                                         income_adequacy.y)),na.rm = TRUE)
describe(ipa_ind$income)
describeBy(ipa_ind$income, ipa_ind$eaipa_country)

cor.test(ipa_ind$income, ipa_ind$age)
cor.test(ipa_ind$income, ipa_ind$degree)

# Security indicator

ipa_ind$security <- rowMeans(subset(ipa_ind, select = c(fin_secu.y, 
                                                       job_stabi.y,
                                                       pension_secu.y)),na.rm = TRUE)
describe(ipa_ind$security)
describeBy(ipa_ind$security, ipa_ind$eaipa_country)

# Correlations with gender, age and education level

cor.test(ipa_ind$fin_secu, ipa_ind$sex)
cor.test(ipa_ind$pension_secu, ipa_ind$age)
cor.test(ipa_ind$pension_secu, ipa_ind$degree)

#exploitation indicator

#revers coloms

cols=c("unpaid_work.y", "pos_toward_overwork.y")
ipa_ind [ ,cols] = 6 - ipa_ind[ ,cols]

ipa_ind$exploit <- rowMeans(subset(ipa_ind, select = c(unpaid_work.y,
                                                       pos_toward_overwork.y)),na.rm = TRUE)
describe(ipa_ind$exploit)
describeBy(ipa_ind$exploit, ipa_ind$eaipa_country)

# Meaning indicator

ipa_ind$meaning <- rowMeans(subset(ipa_ind, select = c(meaningfullness_o_work.y,
                                                       ident_w_prof.y)),na.rm = TRUE)

describe(ipa_ind$meaningfullness_o_work.y)
describeBy(ipa_ind$meaningfullness_o_work.y, ipa_ind$eaipa_country)

describe(ipa_ind$ident_w_prof.y)
describeBy(ipa_ind$ident_w_prof.y, ipa_ind$eaipa_country)

# Correlationtest indicators against each other
cor.test(ipa_ind$income, ipa_ind$security)
cor.test(ipa_ind$security, ipa_ind$exploit)
cor.test(ipa_ind$income, ipa_ind$exploit)
cor.test(ipa_ind$exploit, ipa_ind$meaning)

# CORONA PANDEMIC

# Existential thread to livelihood ----
hist(ipa_ind$unempl_corona)
tabyl(ipa_ind$unempl_corona, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind, unempl_corona, eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind, unempl_corona, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# corona-specific government assistance ---- 

hist(ipa_ind$gov_sup_corona)
tabyl(ipa_ind$gov_sup_corona, ) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind, gov_sup_corona, eu_regions) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

tabyl(ipa_ind, gov_sup_corona, eaipa_country) %>%
  adorn_percentages("col") %>%
  adorn_pct_formatting(digits = 1)

# income differance 2019 / 2020 ----

# Normal distribution test

dfferenz <- ipa_money$income19_eur.x - ipa_money$income20_eur.x
hist(dfferenz)

# T-test for connected samples

t.test(ipa_money$income19_eur.x, ipa_money$income20_eur.x, paired = TRUE, alternative = "greater")
# p value < 0.001 --> Alt hyp true: there is an income decrease of 3150.20 in avrage 

#EU

East <- subset(ipa_money, eu_regions == "Eastern Europe")
dfferenz <- East$income19_eur.x - East$income20_eur.x
hist(dfferenz)
t.test(East$income19_eur.x, East$income20_eur.x, paired = TRUE, alternative = "greater")

West <- subset(ipa_money, eu_regions == "Western Europe")
dfferenz <- West$income19_eur.x - West$income20_eur.x
hist(dfferenz)
t.test(West$income19_eur.x, West$income20_eur.x, paired = TRUE, alternative = "greater")

North <- subset(ipa_money, eu_regions == "Northern Europe")
dfferenz <- North$income19_eur.x - North$income20_eur.x
hist(dfferenz)
t.test(North$income19_eur.x, North$income20_eur.x, paired = TRUE, alternative = "greater")

South <- subset(ipa_money, eu_regions == "Southern Europe")
dfferenz <- South$income19_eur.x - South$income20_eur.x
hist(dfferenz)
t.test(South$income19_eur.x, South$income20_eur.x, paired = TRUE, alternative = "greater")

#Country

Austria <- subset(ipa_money, eaipa_country == "Austria")
dfferenz <- Austria$income19_eur.x - Austria$income20_eur.x
hist(dfferenz)
t.test(Austria$income19_eur.x, Austria$income20_eur.x, paired = TRUE, alternative = "greater")

Germany <- subset(ipa_money, eaipa_country == "Germany")
dfferenz <- Germany$income19_eur.x - Germany$income20_eur.x
hist(dfferenz)
t.test(Germany$income19_eur.x, Germany$income20_eur.x, paired = TRUE, alternative = "greater")

Slovenia <- subset(ipa_money, eaipa_country == "Slovenia")
dfferenz <- Slovenia$income19_eur.x - Slovenia$income20_eur.x
hist(dfferenz)
t.test(Slovenia$income19_eur.x, Slovenia$income20_eur.x, paired = TRUE, alternative = "greater")

Bulgaria <- subset(ipa_money, eaipa_country == "Bulgaria")
dfferenz <- Bulgaria$income19_eur.x - Bulgaria$income20_eur.x
hist(dfferenz)
t.test(Bulgaria$income19_eur.x, Bulgaria$income20_eur.x, paired = TRUE, alternative = "greater") 

Italy <- subset(ipa_money, eaipa_country == "Italy")
dfferenz <- Italy$income19_eur.x - Italy$income20_eur.x
hist(dfferenz)
t.test(Italy$income19_eur.x, Italy$income20_eur.x, paired = TRUE, alternative = "greater") 

Iceland <- subset(ipa_money, eaipa_country == "Iceland")
dfferenz <- Iceland$income19_eur.x - Iceland$income20_eur.x
hist(dfferenz)
t.test(Iceland$income19_eur.x, Iceland$income20_eur.x, paired = TRUE, alternative = "greater")

Sweden <- subset(ipa_money, eaipa_country == "Sweden")
dfferenz <- Sweden$income19_eur.x - Sweden$income20_eur.x
hist(dfferenz)
t.test(Sweden$income19_eur.x, Sweden$income20_eur.x, paired = TRUE, alternative = "greater")

Finland <- subset(ipa_money, eaipa_country == "Finland")
dfferenz <- Finland$income19_eur.x - Finland$income20_eur.x
hist(dfferenz)
t.test(Finland$income19_eur.x, Finland$income20_eur.x, paired = TRUE, alternative = "greater")

Switzerland <- subset(ipa_money, eaipa_country == "Switzerland")
dfferenz <- Switzerland$income19_eur.x - Switzerland$income20_eur.x
hist(dfferenz)
t.test(Switzerland$income19_eur.x, Switzerland$income20_eur.x, paired = TRUE, alternative = "greater")

Romania <- subset(ipa_money, eaipa_country == "Romania")
dfferenz <- Romania$income19_eur.x - Romania$income20_eur.x
hist(dfferenz)
t.test(Romania$income19_eur.x, Romania$income20_eur.x, paired = TRUE, alternative = "greater")

CzechRepublic <- subset(ipa_money, eaipa_country == "Czech Republic")
dfferenz <- CzechRepublic$income19_eur.x - CzechRepublic$income20_eur.x
hist(dfferenz)
t.test(CzechRepublic$income19_eur.x, CzechRepublic$income20_eur.x, paired = TRUE, alternative = "greater")

Hungary <- subset(ipa_money, eaipa_country == "Hungary")
dfferenz <- Hungary$income19_eur.x - Hungary$income20_eur.x
hist(dfferenz)
t.test(Hungary$income19_eur.x, Hungary$income20_eur.x, paired = TRUE, alternative = "greater")

# Decline commissions

ipa_ind <- read_csv("Data/Data_out/ipa_ind.csv")

#2019

describe(ipa_ind$nr_clients19)
prozent <- round(table(ipa_ind$nr_clients19) / sum(table(ipa_ind$nr_clients19)) * 100, 2)
Beschriftung <- c("≤ 3 clients", "3-6 clients", "≥ 7 clients")
pie(table(ipa_ind$nr_clients19), labels = Beschriftung, col=gray(seq(0.3, 1.0, length = 3)))

#2020

describe(ipa_ind$nr_clients20)
prozent <- round(table(ipa_ind$nr_clients20) / sum(table(ipa_ind$nr_clients20)) * 100, 2)
Beschriftung <- c("≤ 3 clients", "3-6 clients", "≥ 7 clients")
pie(table(ipa_ind$nr_clients20), labels = Beschriftung, col=gray(seq(0.3, 1.0, length = 3)))

wilcox.test(ipa_ind$nr_clients19, ipa_ind$nr_clients20, paired = TRUE, exact=FALSE, correct=TRUE, conf.int=TRUE, alternative = "greater")

#EU

describe(East$nr_clients19)
describe(East$nr_clients20)
boxplot(East$nr_clients19, East$nr_clients20)
wilcox.test(East$nr_clients19, East$nr_clients20, paired = TRUE, exact = FALSE, correct = TRUE, conf.int = TRUE, alternative = "greater")

describe(West$nr_clients19)
describe(West$nr_clients20)
boxplot(West$nr_clients19, West$nr_clients20)
wilcox.test(West$nr_clients19, West$nr_clients20, paired = TRUE, exact = FALSE, correct = TRUE, conf.int = TRUE, alternative = "greater")

describe(North$nr_clients19)
describe(North$nr_clients20)
boxplot(North$nr_clients19, North$nr_clients20)
wilcox.test(North$nr_clients19, North$nr_clients20, paired = TRUE, exact = FALSE, correct = TRUE, conf.int = TRUE, alternative = "greater")

describe(South$nr_clients19)
describe(South$nr_clients20)
boxplot(South$nr_clients19, South$nr_clients20)
wilcox.test(South$nr_clients19, South$nr_clients20, paired = TRUE, exact = FALSE, correct = TRUE, conf.int = TRUE, alternative = "greater")

# Decline employment contracts

differenz <- ipa_org$full_time19 - ipa_org$full_time20
hist(differenz)

describe(ipa_org$full_time19)
describe(ipa_org$full_time20)
t.test(ipa_org$full_time19, ipa_org$full_time20, paired = TRUE, alternative = "greater")

  West_org <- subset(ipa_org, eu_regions == "Western Europe")
  describe(West_org$full_time19)
  describe(West_org$full_time20)
  t.test(West_org$full_time19, West_org$full_time20, paired = TRUE, alternative = "greater")

  Eest_org <- subset(ipa_org, eu_regions == "Eastern Europe")
  describe(Eest_org$full_time19)
  describe(Eest_org$full_time20)
  t.test(Eest_org$full_time19, Eest_org$full_time20, paired = TRUE, alternative = "greater")
  
  North_org <- subset(ipa_org, eu_regions == "Northern Europe")
  describe(North_org$full_time19)
  describe(North_org$full_time20)
  t.test(North_org$full_time19, North_org$full_time20, paired = TRUE, alternative = "greater")
  
  South_org <- subset(ipa_org, eu_regions == "Southern Europe")
  describe(South_org$full_time19)
  describe(South_org$full_time20)
  t.test(South_org$full_time19, South_org$full_time20, paired = TRUE, alternative = "greater")

describe(ipa_org$part_time19)
describe(ipa_org$part_time20)
t.test(ipa_org$part_time19, ipa_org$part_time20, paired = TRUE, alternative = "greater")

describe(ipa_org$temporary19)
describe(ipa_org$temporary20)
t.test(ipa_org$temporary19, ipa_org$temporary20, paired = TRUE, alternative = "greater")

describe(ipa_org$fee19)
describe(ipa_org$fee20)
t.test(ipa_org$fee19, ipa_org$fee20, paired = TRUE, alternative = "greater")

  West_org <- subset(ipa_org, eu_regions == "Western Europe")
  describe(West_org$fee19)
  describe(West_org$fee20)
  t.test(West_org$fee19, West_org$fee20, paired = TRUE, alternative = "greater")

  Eest_org <- subset(ipa_org, eu_regions == "Eastern Europe")
  describe(Eest_org$fee19)
  describe(Eest_org$fee20)
  t.test(Eest_org$fee19, Eest_org$fee20, paired = TRUE, alternative = "greater")

  North_org <- subset(ipa_org, eu_regions == "Northern Europe")
  describe(North_org$fee19)
  describe(North_org$fee20)
  t.test(North_org$fee19, North_org$fee20, paired = TRUE, alternative = "greater")

  South_org <- subset(ipa_org, eu_regions == "Southern Europe")
  describe(South_org$fee19)
  describe(South_org$fee20)
  t.test(South_org$fee19, South_org$fee20, paired = TRUE, alternative = "greater")

describe(ipa_org$special_order19)
describe(ipa_org$special_order20)
t.test(ipa_org$special_order19, ipa_org$special_order20, paired = TRUE, alternative = "greater")

# Decline guest performances ----

ipa_guestperf <- ipa_org[ , c(75, 81, 121, 122)]
ipa_guestperf <- subset(ipa_guestperf, guestperf19 < 500)
ipa_guestperf <- subset(ipa_guestperf, guestperf20 < 500)

ipa_guestperf <- na.omit(ipa_guestperf)
dfferenz <- ipa_guestperf$guestperf19 - ipa_guestperf$guestperf20
hist(dfferenz)
describe(ipa_guestperf$guestperf19)
describe(ipa_guestperf$guestperf20)
t.test(ipa_guestperf$guestperf19, ipa_guestperf$guestperf20, paired = TRUE, alternative = "greater")

  West_org <- subset(ipa_guestperf, eu_regions == "Western Europe")
  describe(West_org$guestperf19)
  describe(West_org$guestperf20)
  t.test(West_org$guestperf19, West_org$guestperf20, paired = TRUE, alternative = "greater")

  Eest_org <- subset(ipa_guestperf, eu_regions == "Eastern Europe")
  describe(Eest_org$guestperf19)
  describe(Eest_org$guestperf20)
  t.test(Eest_org$guestperf19, Eest_org$guestperf20, paired = TRUE, alternative = "greater")

  North_org <- subset(ipa_guestperf, eu_regions == "Northern Europe")
  describe(North_org$guestperf19)
  describe(North_org$guestperf20)
  t.test(North_org$guestperf19, North_org$guestperf20, paired = TRUE, alternative = "greater")

  South_org <- subset(ipa_guestperf, eu_regions == "Southern Europe")
  describe(South_org$guestperf19)
  describe(South_org$guestperf20)
  t.test(South_org$guestperf19, South_org$guestperf20, paired = TRUE, alternative = "greater")

# Declines spectators ----

ipa_spec <- ipa_org[ , c(76, 82, 121, 122)]
ipa_spec <- na.omit(ipa_spec)
dfferenz <- ipa_spec$spectators19 - ipa_spec$spectators20
hist(dfferenz)
describe(ipa_spec$spectators19)
describe(ipa_spec$spectators20)
t.test(ipa_spec$spectators19, ipa_spec$spectators20, paired = TRUE, alternative = "greater")

# spectators per region

  West_org <- subset(ipa_spec, eu_regions == "Western Europe")
  describe(West_org$spectators19)
  describe(West_org$spectators20)
  t.test(West_org$spectators19, West_org$spectators20, paired = TRUE, alternative = "greater")
  
  Eest_org <- subset(ipa_spec, eu_regions == "Eastern Europe")
  describe(Eest_org$spectators19)
  describe(Eest_org$spectators20)
  t.test(Eest_org$spectators19, Eest_org$spectators20, paired = TRUE, alternative = "greater")
  
  North_org <- subset(ipa_spec, eu_regions == "Northern Europe")
  describe(North_org$spectators19)
  describe(North_org$spectators20)
  t.test(North_org$spectators19, North_org$spectators20, paired = TRUE, alternative = "greater")
  
  South_org <- subset(ipa_spec, eu_regions == "Southern Europe")
  describe(South_org$spectators19)
  describe(South_org$spectators20)
  t.test(South_org$spectators19, South_org$spectators20, paired = TRUE, alternative = "greater")

# Declines premieres ----
  
ipa_prem <- ipa_org[ , c(71, 77, 121, 122)]  
ipa_prem <- na.omit(ipa_prem)
dfferenz <- ipa_prem$prem19 - ipa_prem$prem20
hist(dfferenz)
describe(ipa_prem$prem19)
describe(ipa_prem$prem20)
t.test(ipa_prem$prem19, ipa_prem$prem20, paired = TRUE, alternative = "greater")

  West <- subset(ipa_prem, eu_regions == "Western Europe")
  describe(West$prem19)
  describe(West$prem20)
  t.test(West$prem19, West$prem20, paired = TRUE, alternative = "greater")

  Eest <- subset(ipa_prem, eu_regions == "Eastern Europe")
  describe(Eest$prem19)
  describe(Eest$prem20)
  t.test(Eest$prem19, Eest$prem20, paired = TRUE, alternative = "greater")

  North <- subset(ipa_prem, eu_regions == "Northern Europe")
  describe(North$prem19)
  describe(North$prem20)
  t.test(North$prem19, North$prem20, paired = TRUE, alternative = "greater")

  South <- subset(ipa_prem, eu_regions == "Southern Europe")
  describe(South$prem19)
  describe(South$prem20)
  t.test(South$prem19, South$prem20, paired = TRUE, alternative = "greater")

# Declines number of  Shows ----

ipa_shows <- ipa_org[ , c(74, 80, 121, 122)]  
ipa_shows <- na.omit(ipa_shows)
dfferenz <- ipa_shows$shows19 - ipa_shows$shows20
hist(dfferenz)
describe(ipa_shows$shows19)
describe(ipa_shows$shows20)
t.test(ipa_shows$shows19, ipa_shows$shows20, paired = TRUE, alternative = "greater")
  
  West <- subset(ipa_shows, eu_regions == "Western Europe")
  describe(West$shows19)
  describe(West$shows20)
  t.test(West$shows19, West$shows20, paired = TRUE, alternative = "greater")
  
  Eest <- subset(ipa_shows, eu_regions == "Eastern Europe")
  describe(Eest$shows19)
  describe(Eest$shows20)
  t.test(Eest$shows19, Eest$shows20, paired = TRUE, alternative = "greater")
  
  North <- subset(ipa_shows, eu_regions == "Northern Europe")
  describe(North$shows19)
  describe(North$shows20)
  t.test(North$shows19, North$shows20, paired = TRUE, alternative = "greater")
  
  South <- subset(ipa_shows, eu_regions == "Southern Europe")
  describe(South$shows19)
  describe(South$shows20)
  t.test(South$shows19, South$shows20, paired = TRUE, alternative = "greater")
  
# Declines shows in avr ----

ipa_avr <- ipa_org[ , c(73, 79, 121, 122)]  
ipa_avr <- na.omit(ipa_avr)
dfferenz <- ipa_avr$shows_avre19 - ipa_avr$shows_avr20
hist(dfferenz)
describe(ipa_avr$shows_avre19)
describe(ipa_avr$shows_avr20)
t.test(ipa_avr$shows_avre19, ipa_avr$shows_avr20, paired = TRUE, alternative = "greater")
  
  West <- subset(ipa_avr, eu_regions == "Western Europe")
  describe(West$shows_avre19)
  describe(West$shows_avr20)
  t.test(West$shows_avre19, West$shows_avr20, paired = TRUE, alternative = "greater")
  
  Eest <- subset(ipa_avr, eu_regions == "Eastern Europe")
  describe(Eest$shows_avre19)
  describe(Eest$shows_avr20)
  t.test(Eest$shows_avre19, Eest$shows_avr20, paired = TRUE, alternative = "greater")
  
  North <- subset(ipa_avr, eu_regions == "Northern Europe")
  describe(North$shows_avre19)
  describe(North$shows_avr20)
  t.test(North$shows_avre19, North$shows_avr20, paired = TRUE, alternative = "greater")
  
  South <- subset(ipa_avr, eu_regions == "Southern Europe")
  describe(South$shows_avre19)
  describe(South$shows_avr20)
  t.test(South$shows_avre19, South$shows_avr20, paired = TRUE, alternative = "greater")
    
  # Declines productions ----

  prod <- ipa_org[ , c(72, 78, 121, 122)]  
  prod <- na.omit(prod)
  dfferenz <- prod$prod19 - prod$prod20
  hist(dfferenz)
  describe(prod$prod19)
  describe(prod$prod20)
  t.test(prod$prod19, prod$prod20, paired = TRUE, alternative = "greater")
  
  West <- subset(prod, eu_regions == "Western Europe")
  describe(West$prod19)
  describe(West$prod20)
  t.test(West$prod19, West$prod20, paired = TRUE, alternative = "greater")
  
  Eest <- subset(prod, eu_regions == "Eastern Europe")
  describe(Eest$prod19)
  describe(Eest$prod20)
  t.test(Eest$prod19, Eest$prod20, paired = TRUE, alternative = "greater")
  
  North <- subset(prod, eu_regions == "Northern Europe")
  describe(North$prod19)
  describe(North$prod20)
  t.test(North$prod19, North$prod20, paired = TRUE, alternative = "greater")
  
  South <- subset(prod, eu_regions == "Southern Europe")
  describe(South$prod19)
  describe(South$prod20)
  t.test(South$prod19, South$prod20, paired = TRUE, alternative = "greater")

#GOVERNANCE

library(MVN)
library(lavaan)
library(semPlot)

ipa_common <- read_csv("Data/Data_out/ipa_common.csv")
ipa_ind <- read_csv("Data/Data_out/ipa_ind.csv")
ipa_org <- read_csv("Data/Data_out/ipa_org.csv")
ipa_genre_all <- read_csv("Data/Data_out/_genre/ipa_genre_all.csv")
ipa_profession <- read_csv("Data/Data_out/_prof_org/ipa_profession.csv")
ipa_wage <- read_csv("Data/Data_out/_money/ipa_wage.csv")
exp_data <- read_excel("~/Documents/_PhD Promotion/Analyses/Qual Analyse/Daten/data_Interview_IPA_Association_2021-04-19_13-25.xlsx")

# Membership EAIPA member organization

tabyl(ipa_common, member_eaipa, )
tabyl(ipa_common, member, ) 

Austria <- subset(ipa_common, eaipa_country == "Austria")
tabyl(Austria$member_eaipa, ) 
tabyl(Austria$member_other, ) 

Bulgaria <- subset(ipa_common, eaipa_country == "Bulgaria")
tabyl(Bulgaria$member_eaipa, ) 
tabyl(Bulgaria$member_other, ) 

CzechRepublic <- subset(ipa_common, eaipa_country == "Czech Republic")
tabyl(CzechRepublic$member_eaipa, )
tabyl(CzechRepublic$member_other, ) 

Finland <- subset(ipa_common, eaipa_country == "Finland")
tabyl(Finland$member_eaipa, )
tabyl(Finland$member_other, ) 

Germany <- subset(ipa_common, eaipa_country == "Germany")
tabyl(Germany$member_eaipa, )
tabyl(Germany$member_other, ) 

Hungary <- subset(ipa_common, eaipa_country == "Hungary")
tabyl(Hungary$member_eaipa, )
tabyl(Hungary$member_other, ) 

Iceland <- subset(ipa_common, eaipa_country == "Iceland")
tabyl(Iceland$member_eaipa, )
tabyl(Iceland$member_other, ) 

Italy <- subset(ipa_common, eaipa_country == "Italy")
tabyl(Italy$member_eaipa, )
tabyl(Italy$member_other, ) 

Romania <- subset(ipa_common, eaipa_country == "Romania")
tabyl(Romania$member_eaipa, )
tabyl(Romania$member_other, ) 

Slovenia <- subset(ipa_common, eaipa_country == "Slovenia")
tabyl(Slovenia$member_eaipa, )
tabyl(Slovenia$member_other, ) 

Sweden <- subset(ipa_common, eaipa_country == "Sweden")
tabyl(Sweden$member_eaipa, )
tabyl(Sweden$member_other, ) 

Switzerland <- subset(ipa_common, eaipa_country == "Switzerland")
tabyl(Switzerland$member_eaipa, )
tabyl(Bulgaria$member_other, ) 

# Representation and democratic organization ----
# Filtering members only

ipa_member <- subset(ipa_common, member_eaipa > 0 )

# Representation 

describe(ipa_member$rep_through_advo)
describeBy(ipa_member$rep_through_advo, ipa_member$eu_regions)
describeBy(ipa_member$rep_through_advo, ipa_member$eaipa_country)

# Democratic organization

describe(ipa_member$demo_asso)
describeBy(ipa_member$demo_asso, ipa_member$eu_regions)
describeBy(ipa_member$demo_asso, ipa_member$eaipa_country)

# Active participation 

describe(ipa_member$particip_i_advo)
describeBy(ipa_member$particip_i_advo, ipa_member$eu_regions)
describeBy(ipa_member$particip_i_advo, ipa_member$eaipa_country)

#Build governance scale 

ipa_gov <- subset(ipa_common, select = c(case, policy_feedback_pos:dialog_gov_nat, coop_gov, policy_inclu))

str(ipa_gov)
summary(ipa_gov)

# % NA in row

percentmissing = function(x){ sum(is.na(x))/length(x)*100}
missing = apply(ipa_gov, 1, percentmissing)
table(missing)

Delete rows with more the 20% NA

replaceorg = subset(ipa_gov, missing <= 20)

# % NA in column 

apply(replaceorg, 2, percentmissing)

# Multiple impute datapoints

tempnomiss = mice(replaceorg) 
nomiss = complete(tempnomiss, 1)
summary(nomiss)

# Build scale 

governance = c(2, 3, 4, 5, 6, 7) 

# Test reliability 

psych::alpha(ipa_gov[governance])

ipa_gov$governance = apply(ipa_gov[ ,governance], 1, mean)

# Merge Datasets

ipa_common_final <- merge( ipa_common, ipa_gov, by = "case")

# Describe governance index

describe(ipa_common_final$governance)

# Boxplot governance per region / country

boxplot(ipa_common_final$governance ~ ipa_common_final$eu_regions)

ggplot(data = ipa_common_final, mapping = aes(x=governance, 
                                       y=eu_regions)) + 
  geom_boxplot() +
  theme_linedraw() + 
  theme(legend.position = "bottom") + 
  scale_fill_brewer(palette="Greys") +
  labs(title = "",
       x = "Level of perceived governance opportunity",
       y = "",
       fill = "",
       caption = "") 

ggplot(data = ipa_common_final, mapping = aes(x=governance, 
                                              y=eaipa_country)) + 
  geom_boxplot() +
  theme_linedraw() + 
  theme(legend.position = "bottom") + 
  scale_fill_brewer(palette="Greys") +
  labs(title = "",
       x = "Level of perceived governance opportunity",
       y = "",
       fill = "",
       caption = "") 

# Describe governance per region / country

describeBy(ipa_common_final$governance, ipa_common_final$eu_regions)
describeBy(ipa_common_final$governance, ipa_common_final$eaipa_country)

# Information by government

describe(ipa_common$policy_com_gov)
describeBy(ipa_common$policy_com_gov, ipa_common$eu_regions)
describeBy(ipa_common$policy_com_gov, ipa_common$eaipa_country)

# Consultation and dialogue

ipa_common$dialog <- rowMeans(subset(ipa_common, select = c(policy_feedback_pos, dialog_gov_mun, dialog_gov_nat)),na.rm = TRUE)
describe(ipa_common$dialog)
describeBy(ipa_common$dialog, ipa_common$eu_regions)
describeBy(ipa_common$dialog, ipa_common$eaipa_country)

# Codetermination

ipa_common$coditermination <- rowMeans(subset(ipa_common, select = c(policy_feedback_consid, coop_gov, policy_inclu)),na.rm = TRUE)
describe(ipa_common$coditermination)
describeBy(ipa_common$coditermination, ipa_common$eu_regions)
describeBy(ipa_common$coditermination, ipa_common$eaipa_country)

# Assessment of experts

# Information by government

exp_data$information <- rowMeans(subset(exp_data, select = c(pa03_01)),na.rm = TRUE)
describe(exp_data$information)
describeBy(exp_data$information, exp_data$eu_regions)
describeBy(exp_data$information, exp_data$a002_02)

# Consultation and dialogue

exp_data$dialog <- rowMeans(subset(exp_data, select = c(pa03_03, pa03_05)),na.rm = TRUE)
describe(exp_data$dialog)
describeBy(exp_data$dialog, exp_data$eu_regions)
describeBy(exp_data$dialog, exp_data$a002_02)

# Codetermination

exp_data$coditermination <- rowMeans(subset(exp_data, select = c(pa03_06)),na.rm = TRUE)
describe(exp_data$coditermination)
describeBy(exp_data$coditermination, exp_data$eu_regions)
describeBy(exp_data$coditermination, exp_data$a002_02)

# Create governance indicator

exp_data$governance <- rowMeans(subset(exp_data, select = c(pa03_01, 
                                                                pa03_02,
                                                                pa03_03,
                                                                pa03_05, 
                                                                pa03_06)),na.rm = TRUE)

# Plot Expert Governance assessment

ggplot(ipa_common, aes(x=representation, y=eaipa_country)) + 
  geom_boxplot()

ggplot(ipa_common, aes(x=representation, y=eaipa_country)) + 
  geom_boxplot()

# Cultural policy issues of particular relevance (Data from expert Questionnaire) ----

exp_data <- read_excel("~/Documents/_PhD Promotion/Analyses/Qual Analyse/Daten/data_Interview_IPA_Association_2021-04-19_13-25.xlsx")

# Clean data

exp_data <- rename(exp_data,
                   proj_fun_nat = RA06_01,
                   schol_nat = RA06_02,
                   struct_fun_nat = RA06_03,
                   newc_fun_nat = RA06_04, 
                   mobi_nat = RA06_05,
                   int_mobi_nat = RA06_06,
                   venue_fund_nat = RA06_07,
                   proj_fun_mun = RA07_01,
                   schol_mun = RA07_02,
                   struct_fun_mun = RA07_03,
                   newc_fun_mun = RA07_04, 
                   mobi_mun = RA07_05,
                   int_mobi_mun = RA07_06,
                   venue_fund_mun = RA07_07,
                   
                   focus_pol = PA01_01,
                   focus_cons = PA01_02,
                   focus_net = PA01_03,
                   focus_prod = PA01_04, 
                   
                   interest_memberdialog = PA02_01,
                   interest_memberneeds = PA02_02,
                   interest_imp_fund_struc = PA02_03,
                   interest_incr_fund = PA02_04,
                   interest_imp_work_cond = PA02_05,
                   interest_create_qual_offers = PA02_06,
                   interest_create_pres_opp = PA02_07,
                   interest_create_distr_opp = PA02_08,
                   interest_create_netw_opp = PA02_09)

exp_data <- clean_names(exp_data)

exp_data <- exp_data[-c(8,9,13), ]

exp_data$eu_regions[exp_data$a002_02 %in% c("Austria",
                                            "Germany",
                                            "Switzerland")] <- "Western Europe"

exp_data$eu_regions[exp_data$a002_02 %in% c("Bulgaria",
                                            "Czech Republic",
                                            "Hungary",
                                            "Romania")] <- "Eastern Europe"

exp_data$eu_regions[exp_data$a002_02 %in% c("Italy",
                                            "Slovenia")] <- "Southwestern Europe"

exp_data$eu_regions[exp_data$a002_02 %in% c("Finland",
                                            "Iceland",
                                            "Sweden")] <- "Scandinavia"

# Describe interests

describe(exp_data$interest_incr_fund)
describe(exp_data$interest_imp_fund_struc)
describe(exp_data$interest_imp_work_cond)

describe(exp_data$interest_memberneeds)
describe(exp_data$interest_create_netw_opp)
describe(exp_data$interest_memberdialog)
describe(exp_data$interest_create_distr_opp)
describe(exp_data$interest_create_qual_offers)
describe(exp_data$interest_create_pres_opp)