Difference between revisions of "Goherr: Fish consumption study"

  
library(OpasnetUtils)
library(ggplot2)
library(reshape2)

objects.latest("Op_en7749", "preprocess")

for(i in c(5:6, 16, 29:30, 46:49, 85:86, 95:98, 135)) {
  temp <- survey[!is.na(survey[[i]]),]
    p <- ggplot(temp, aes(x = temp[[i]])) + geom_bar() +
    theme_gray(base_size = 18) + theme(axis.text.x=element_text(angle=90,hjust=1,vjust=0.4)) +
    labs(title = colnames(temp[i])) + xlab("") + facet_wrap(~ Country)
  print(p)
}

temp <- melt(survey, measure.vars = 145:153, variable.name = "Changing.factor", value.name = "Impact")
levs <- c("-5 strongly disagree", "-4", "-3", "-2", "-1", "0 Neutral", "1", "2", "3", "4", "5 strongly agree", "I don't know")
temp$Impact <- factor(temp$Impact, levels = levs, labels = levs, ordered = TRUE)

ggplot(temp, aes(x = Impact, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Changing.factor)+
  theme_gray(base_size = 24) + theme(axis.text.x=element_text(angle=90,hjust=1,vjust=0.4))

temp <- melt(survey, measure.vars = 99:108, variable.name = "Reason.to.eat", value.name = "Impact")
ggplot(temp, aes(x = Impact, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Reason.to.eat)+
   theme_gray(base_size = 24)

temp <- melt(survey, measure.vars = 38:43, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Baltic salmon sources")

temp <- melt(survey, measure.vars = 50:59, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Reasons to eat Baltic salmon")

temp <- melt(survey, measure.vars = 62:70, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Reasons not to eat Baltic salmon")

temp <- melt(survey, measure.vars = 73:82, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Baltic salmon actions")

temp <- melt(survey, measure.vars = 87:92, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Baltic herring sources")

temp <- melt(survey, measure.vars = 99:108, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Reasons to eat Baltic herring")

temp <- melt(survey, measure.vars = 111:120, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Reasons not to eat Baltic herring")

temp <- melt(survey, measure.vars = 123:132, variable.name = "Variable", value.name = "Value")
ggplot(temp, aes(x = Value, fill = Country)) + geom_histogram(position = "dodge") + facet_wrap(~ Variable)+
   theme_gray(base_size = 24) + labs(title = "Baltic herring actions")

library(OpasnetUtils)
library(reshape2)
library(ggplot2)

############## Generic functions and objects are defined first.

### webropol.convert converts a csv file from Webropol into a useful data.frame.

webropol.convert <- function(
  data, # Data.frame created from a Webropol csv file. The first row should contain headings.
  rowfact, # Row number where the factor levels start (in practice, last row + 3)
  textmark = "Other open" # The text that is shown in the heading if there is an open sub-question.
) {
  out <- dropall(data[2:(rowfact - 3) , ])
  subquestion <- t(data[1 , ])
  subquestion <- gsub("\xa0", " ", subquestion)
  subquestion <- gsub("\xb4", " ", subquestion)
  subquestion <- gsub("\n", " ", subquestion)
  #  subquestion <- gsub("\\(", " ", subquestion)
  #  subquestion <- gsub("\\)", " ", subquestion)
  textfield <- regexpr(textmark, subquestion) != -1
  subquestion <- strsplit(subquestion, ":") # Divide the heading into a main question and a subquestion.
  subqtest <- 0 # The previous question name.
  for(i in 1:ncol(out)) {
    #print(i)
    if(subquestion[[i]][1] != subqtest) { # If part of previous question, use previous fact.
      fact <- as.character(data[rowfact:nrow(data) , i]) # Create factor levels from the end of Webropol file.
      fact <- fact[fact != ""] # Remove empty rows
      fact <- gsub("\xa0", " ", fact)
      fact <- gsub("\xb4", " ", fact)
      fact <- gsub("\n", " ", fact)
      fact <- strsplit(fact, " = ") # Separate value (level) and interpretation (label)
    }
    if(length(fact) != 0 & !textfield[i]) { # Do this only if the column is not a text type column.
      out[[i]] <- factor(
        out[[i]], 
        levels = unlist(lapply(fact, function(x) x[1])), 
        labels = unlist(lapply(fact, function(x) x[2])), 
        ordered = TRUE
      )
    }
    subqtest <- subquestion[[i]][1]
  }
  return(out)
}

# merge.questions takes a multiple checkbox question and merges that into a single factor.
# First levels in levs have priority over others, if several levels apply to a row.

merge.questions <- function(
  dat, # data.frame with questionnaire data
  cols, # list of vectors of column names or numbers to be merged into one level in the factor
  levs, # vector (with the same length as cols) of levels of factors into which questions are merged.
  name # text string for the name of the new factor column in the data.
) {
  for(i in length(cols):1) {
    temp <- FALSE
    for(j in rev(cols[[i]])) {
      temp <- temp | !is.na(dat[[j]])
    }
    dat[[name]][temp] <- levs[i]
  }
  dat[[name]] <- factor(dat[[name]], levels = levs, ordered = TRUE)
  return(dat)
}

############# Data preprocessing

# Get the data either from Opasnet or your own hard drive.

#Survey original file: N:/Ymal/Projects/Goherr/WP5/Goherr_fish_consumption.csv

survey <- opasnet.csv("5/57/Goherr_fish_consumption.csv", sep = ";", fill = TRUE, quote = "\"")
surcol <- t(survey[1,])

survey <- webropol.convert(survey, 2121, textmark = ":Other open") # Data file is converted to data.frame using levels at row 1269.

# Take the relevant columnames from the table on the page.
colframe <- opbase.data("Op_en7749", subset = "Questions in the Goherr questionnaire")
colnames(survey) <- colframe$Result

objects.store(survey, surcol)
cat("objects survey, surcol were stored.\n")

# This is code Op_en7749/ on page [[Goherr: Fish consumption study]]

library(OpasnetUtils)
library(ggplot2)
library(reshape2)

objects.latest("Op_en7749", "preprocess")
assumptions <- opbase.data("Op_en7749.assumptions_for_calculations")

survey$Row <- 1:nrow(survey)
survey$Ages <- ifelse(as.numeric(as.character(survey$Age)) < 46, "18-45",">45")
short <- survey[c("Country", "Gender", "Row", "Ages")]

salmon_dishtimes <- Ovariable(
  "salmon_dishtimes", 
  data = data.frame(
    short,
    Value = survey[[46]],
    Result = 1
  )
) #salmon with dish times / year

salmon_dishamount <- Ovariable(
  "salmon_dishamount", 
  data = data.frame(
    short,
    Value = survey[[47]],
    Result = 1
  )
) #salmon with dish amount per one serving


salmon_dishtimesassump <- Ovariable(
  "salmon_dishtimesassump", 
  data = assumptions[assumptions$Variable == "V46" , c("Value", "Result")]
)

salmon_dishamountassump <- Ovariable(
  "salmon_dishamountassump", 
  data = assumptions[assumptions$Variable == "V47" , c("Value", "Result")]
)

amount <- Ovariable("amount", 
  dependencies = data.frame(Name = c(
    "salmon_dishamount",
    "salmon_dishamountassump",
    "salmon_dishtimes",
    "salmon_dishtimesassump"
  )),
  formula = function(...) {
    
    out <- ((salmon_dishamount * salmon_dishamountassump) * (salmon_dishtimes * salmon_dishtimesassump))/ 365
    # Per vuosi -> per d
    
    #result(out)[result(out) == 0] <- 0.01 # Ei jätetä nollia saantiin
    return(out)
  }
)

temp1 <- EvalOutput(amount, forceEval = TRUE, N = 3)
oprint(head(temp1@output))