This is something I keep looking up, because for whatever reason R does not come with a built-in function to calculate the mode. (The mode() function does something else, not what I’d expect given that there are mean() and median()…) It’s quite easy to write a short function to calculate the mode in R:
Mode <- function(x) {
uni <- unique(x)
uni[which.max(tabulate(match(x, uni)))]
}
Didier Ruedin 
Mode <- function(x) {
uni <- unique(x)
uni[which.max(tabulate(match(x, uni)))]
}
When there is no mode, this function returns a mode e.g.
fruit = c(rep('apple', 5), rep('pear', 5), rep('banana', 5))
Mode(fruit)
I think it is more common to say that there is more than one mode in this case. So yes, the function only gives the first mode if there are multiple ones. I have something for this: library(agrmt); modes(collapse(fruit)) # this is a function I have written to identify multiple modes. The handling is a bit different, because it uses frequency vectors — hence the collapse() function in the middle. Perhaps this helps?