R/AAstat.R

In seqinr: Biological Sequences Retrieval and Analysis

utils::data(sysdata, envir=environment())
AAstat <- function(seq, plot = TRUE){
    #
    # seq is a protein sequence as a vector of (upper case) chars.
    #
    AAP <- SEQINR.UTIL$AA.PROPERTY
    tutu <- lapply(names(AAP), function(x) which(seq %in% AAP[[x]]))
    names(tutu) <- names(AAP)
    
    n.items <- length(tutu) # Number of physoco-chemical properties
    n.res <- length(seq)    # Number of residues in the protein
    
    if(plot == TRUE){
        coul <- grDevices::rainbow(n.items)
        graphics::plot(c(0, n.res), c(0, n.items + 1), type = "n", axes = FALSE, 
             ann = FALSE, xlim = c(0, n.res + 1))
        graphics::title(xlab = "Position of the residues along the sequence")
        graphics::axis(2, at = seq(1.5, 10, 1), labels = names(tutu), col.lab = "blue", las = 1, cex.axis = 0.8)
        graphics::axis(1, at = seq(0, n.res, 15), labels = seq(0, n.res, 15), col.axis = "blue")
        lapply(seq_len(n.items), function(x){
            graphics::segments(tutu[[x]], x, tutu[[x]], x + 1, col = coul[x], lwd = 2)
            graphics::rect(0, x, n.res, 1, lwd = 2)
        })
        graphics::rect(0, n.items + 1, n.res, 1, lwd = 2)
    }
    
    res1 <- lapply(tutu,function(x){length(x)/n.res})
    res2 <- table(factor(seq, levels = levels(SEQINR.UTIL$CODON.AA$L)))
    res3 <- computePI(seq)
    return(list(Compo = res2, Prop = res1, Pi = res3))
}