# HAP40

library(readr)

####################################################################
## Data import from ConSurf                                  #######
####################################################################

data_con <-
  read_delim(
    "Consurf AA Conservation Scores modified.csv",
    ";",
    escape_double = FALSE,
    trim_ws = TRUE
  )


####################################################################
## Caculation of conservation score for each domain             ####
####################################################################

#' Function to average conservation scorce for certain protein region
#' Scores obtained by ConSturf method)
#'
#' @param data: R dataframe containing postions and scores
#' @param pos: vector with start and end position of protein region
#'
#' @return: vector with elements pos, mean, SEM

conservation <- function(data, pos) {
  data <- as.data.frame(data)                                    # Loading of conservation scores
  data.sub <- as.vector(data[pos, 3])                            # Subsetting of data frame
  
  mean.conservation <- mean(data.sub)                            # calculation of mean
  SEM.conservation <- sd(data.sub)/sqrt(length(pos))             # calculation of SEM
  
  # Defining output of function
  output <-c( 
    mean.conservation, 
    SEM.conservation)
  names(output) <- c("mean", "SEM")
  return(output)
}

#####################################################################
## List defining start and end point of the corresponding domains ###
## Note: in this case: notation in Guo et al. 2018 + 6            ###
#####################################################################

# Positions

regions <- list(1:216,217:258,259:371)

## Name of regions
names(regions) <-
  c(
    "N-terminal domain",
    "central",
    "C-terminal domain"
  )

###########################################################################
## Repetition of function conservation for each element of list regions ###
###########################################################################

x <- sapply(regions,FUN = conservation,data=data_con)

library(readr)

####################################################################
## Data import from ConSurf                                  #######
####################################################################

data_con <-
  read_delim(
    "Consurf AA Conservation Scores modified.csv",
    ";",
    escape_double = FALSE,
    trim_ws = TRUE
  )


####################################################################
## Caculation of conservation score for each domain             ####
####################################################################

#' Function to average conservation scorce for certain protein region
#' Scores obtained by ConSturf method)
#'
#' @param data: R dataframe containing postions and scores
#' @param pos: vector with start and end position of protein region
#'
#' @return: vector with elements pos, mean, SEM

conservation <- function(data, pos) {
  data <- as.data.frame(data)                                    # Loading of conservation scores
  data.sub <- as.vector(data[pos, 2])                            # Subsetting of data frame
  
  mean.conservation <- mean(data.sub)                            # calculation of mean
  SEM.conservation <- sd(data.sub)/sqrt(length(pos))             # calculation of SEM
  
  # Defining output of function
  output <-c( 
    mean.conservation, 
    SEM.conservation)
  names(output) <- c("mean", "SEM")
  return(output)
}

#####################################################################
## List defining start and end point of the corresponding domains ###
## Note: in this case: notation in Guo et al. 2018 + 6            ###
#####################################################################

# Positions
C_HEAT.1 <- c(2098:2126,2463:2515,2670:3104)   # Positions of C-HEAT without insertions
N_HEAT.1 <- c(97:405,691:1690)
regions <- list(97:1690,1690:2097,2098:3104,2127:2462,2516:2669,C_HEAT.1,406:680,N_HEAT.1)

## Name of regions
names(regions) <-
  c(
    "N-HEAT",
    "bridge",
    "C-HEAT",
    "C-Insertion_1",
    "C-Insertion_2",
    "C-HEAT without Insertions",
    "N-Insertion",
    "N-HEAT without Insertion"
  )

###########################################################################
## Repetition of function conservation for each element of list regions ###
###########################################################################

x <- sapply(regions,FUN = conservation,data=data_con)
