Resolved: Summarise Data using group_by across several columns

set.seed(123)
gender <- sample(1:2,12,replace=T)
country <- c('FIN', 'FIN', 'EST', 'NIG','NIG','JAM', 'FIN', 'NIG', 'EST', 'NIG','NIG','JAM','FIN', 'FIN', 'EST', 'NIG','NIG','JAM', 'FIN', 'NIG', 'EST', 'NIG','NIG','JAM')
a <- sample(1:5,24,,replace=T)
b <- sample(1:5,24,,replace=T)
c <- sample(1:5,24,,replace=T)
d <- sample(1:5,24,,replace=T)
# Join the variables to create a data frame
df <- data.frame(gender,country,a,b,c,d)

df <- df %>% mutate_at(c("gender","country","a","b","c","d"), as.factor)

df %>% 
  group_by(gender,country,a) %>%
  summarise(totaln=n()) %>% 
  group_by(country,gender) %>% 
  mutate(percentage=totaln/sum(totaln)*100)

   gender country a     totaln percentage
   <fct>  <fct>   <fct>  <int>      <dbl>
 1 1      EST     1          2       50  
 2 1      EST     3          1       25  
 3 1      EST     4          1       25  
 4 1      FIN     1          1       25  
 5 1      FIN     2          3       75  
 6 1      NIG     1          1       25  
 7 1      NIG     2          1       25  
 8 1      NIG     3          1       25  
 9 1      NIG     4          1       25  
10 2      FIN     1          1       50  
11 2      FIN     3          1       50  
12 2      JAM     1          2       50  
13 2      JAM     3          2       50  
14 2      NIG     3          1       16.7
15 2      NIG     4          1       16.7
16 2      NIG     5          4       66.7

df_summary <- df %>% pivot_longer(c(a, b, c, d), names_to = "abcd") %>%
  group_by(gender, country, abcd, value) %>%
  summarise(totaln = n()) %>%
  group_by(country, gender, abcd) %>%
  mutate(percentag = totaln / sum(totaln) * 100)

df_summary 

#> # A tibble: 64 × 6
#> # Groups:   country, gender, abcd [24]
#>    gender country abcd  value totaln percentag
#>    <fct>  <fct>   <chr> <fct>  <int>     <dbl>
#>  1 1      EST     a     1          1        25
#>  2 1      EST     a     3          2        50
#>  3 1      EST     a     4          1        25
#>  4 1      EST     b     2          1        25
#>  5 1      EST     b     3          1        25
#>  6 1      EST     b     4          2        50
#>  7 1      EST     c     2          1        25
#>  8 1      EST     c     5          3        75
#>  9 1      EST     d     1          1        25
#> 10 1      EST     d     2          2        50
#> # … with 54 more rows

df_summary %>% 
  ungroup() %>% 
  mutate(summary = paste0(round(percentag), "% (n = ", totaln, ")")) %>% 
  select(gender, country, abcd, value, summary) %>% 
  pivot_wider(names_from = c(gender, country), values_from = summary, values_fill = "0")