diff --git a/.Rhistory b/.Rhistory
deleted file mode 100644
index 552cf98..0000000
--- a/.Rhistory
+++ /dev/null
@@ -1,512 +0,0 @@
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 45.5, sigmax = 0.01, sigmay = 0.01, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -4, muy = 5.5, sigmax = 0.01, sigmay = 0.01, rhoxy = -0.5))
-max(3,5,2)
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 45.5, sigmax = 7, sigmay = 9, rhoxy = -0.5))
-library(eyeball)
-library(eyeball)
-devtools::load_all(".")
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 45.5, sigmax = 7, sigmay = 9, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 3, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 2.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-library(eyeball)
-library(eyeball)
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-devtools::load_all(".")
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 1, distribution.pairs = "Normal", mux = -1, muy = 4.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.1, distribution.pairs = "Normal", mux = -1, muy = 4.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.1, distribution.pairs = "Normal", mux = -1, muy = 2.5, sigmax = 5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.1, distribution.pairs = "Normal", mux = -1, muy = 2.5, sigmax = 45, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.1, distribution.pairs = "Normal", mux = -1, muy = 2.5, sigmax = 1.5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.1, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 1.5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.1, distribution.pairs = "Normal", mux = -1, muy = 0.25, sigmax = 1.5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.4, distribution.pairs = "Normal", mux = -1, muy = 0.25, sigmax = 1.5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.4, distribution.pairs = "Normal", mux = -1, muy = 0.25, sigmax = 4.5, sigmay = 2, rhoxy = -0.5))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade",S = 500, C= 0.4, distribution.pairs = "Normal", mux = -1, muy = 0.25, sigmax = 4.5, sigmay = 2, rhoxy = -0.5))
-library(eyeball)
-devtools::load_all(".")
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Niche", S = 500, C = 0.25, distribution.pairs = "FourCorner", mux = -5, muy = 2, sigmax = 3, sigmay = 2, rhoxy = -2/3))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Niche", S = 500, C = 0.25, distribution.pairs = "FourCorner", mux = -5, muy = 2, sigmax = 3, sigmay = 2, rhoxy = -2/3))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Niche", S = 500, C = 0.25, distribution.pairs = "FourCorners", mux = -5, muy = 2, sigmax = 3, sigmay = 2, rhoxy = -2/3))
-library(eyeball)
-devtools::load_all(".")
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Niche", S = 500, C = 0.25, distribution.pairs = "FourCorner", mux = -5, muy = 2, sigmax = 3, sigmay = 2, rhoxy = -2/3))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 500, C = 0.25, distribution.pairs = "FourCorner", mux = -5, muy = 2, sigmax = 3, sigmay = 2, rhoxy = -2/3))
-eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 500, C = 0.25, distribution.pairs = "FourCorner", mux = -1, muy = 2, sigmax = 3, sigmay = 2, rhoxy = -2/3))
-library(eyeball)
-devtools::load_all(".")
-library(eyeball)
-devtools::load_all(".")
-Approx <- eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.25, foodweb.model = "Cascade"))
-Approx
-if (is.null(Approx$M.eigenvalues)){
-Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
-}
-## Data frame for the eigenvalues of M
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-eye.ellipse.df <- function(centerx, radiusx, radiusy){
-thetas <- seq(pi / 2.0, 0.0, length = 1000)
-xbase <- radiusx * cos(thetas)
-ybase <- radiusy * sin(thetas)
-x <- c(xbase, rev(xbase), -xbase, rev(-xbase))
-y <- c(ybase, rev(-ybase), -ybase, rev(ybase))
-return(data.frame(Real = x + centerx, Imaginary = y))
-}
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-## Ellipse and extra eigenvalue for Tang et al. criteria
-st <- Approx$TangEtAl.stats
-radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
-radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius.h, radius.v)
-tmp$Type <- "Tang et al."
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "Tang et al."))
-pl
-pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA)
-pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA)
-st <- Approx$eyeball.stats
-st
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(center.A, radius.A, radius.A)
-tmp$Type <- "eyeball, A"
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(center.A, radius.A, radius.A)
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(st$center.A, st$radius.A, st$radius.A)
-tmp$Type <- "eyeball, A"
-ev.polygon <- rbind(ev.polygon, tmp)
-tmp <- eye.ellipse.df(st$radius.A + st$center.A, st$radius.B.h, st$radius.B.v)
-tmp$Type <- "eyeball, B"
-ev.polygon <- rbind(ev.polygon, tmp)
-pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA)
-pl
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 3)
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 4)
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 2)
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1)
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1)
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-pl + geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-only.bulk <- FALSE
-if (is.null(Approx$M.eigenvalues)){
-Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
-}
-## Data frame for the eigenvalues of M
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-## Ellipse and extra eigenvalue for Tang et al. criteria
-st <- Approx$TangEtAl.stats
-radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
-radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius.h, radius.v)
-tmp$Type <- "Tang et al."
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "Tang et al."))
-## Circle and Ellipse for eyeball approximation
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(st$center.A, st$radius.A, st$radius.A)
-tmp$Type <- "eyeball, A"
-ev.polygon <- rbind(ev.polygon, tmp)
-tmp <- eye.ellipse.df(st$radius.A + st$center.A, st$radius.B.h, st$radius.B.v)
-tmp$Type <- "eyeball, B"
-ev.polygon <- rbind(ev.polygon, tmp)
-## Now draw the eigenvalues and the approximations
-pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA) +
-geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-## and cut the plot
-if (only.bulk == FALSE) {
-pl <- pl + coord_cartesian(xlim = 1.1 * range(Re(Approx$M.eigenvalues)),
-ylim = 1.1 * range(Im(Approx$M.eigenvalues)))
-}
-pl
-## Check if the eigenvalues of M have been computed already
-if (is.null(Approx$M.eigenvalues)){
-Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
-}
-## Data frame for the eigenvalues of M
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-## Ellipse and extra eigenvalue for Tang et al. criteria
-st <- Approx$TangEtAl.stats
-radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
-radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius.h, radius.v)
-tmp$Type <- "Tang et al."
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "Tang et al."))
-## Circle and Ellipse for eyeball approximation
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(st$center.A, st$radius.A, st$radius.A)
-tmp$Type <- "eyeball, A"
-ev.polygon <- rbind(ev.polygon, tmp)
-tmp <- eye.ellipse.df(st$radius.A + st$center.A, st$radius.B.h, st$radius.B.v)
-tmp$Type <- "eyeball, shifted B"
-ev.polygon <- rbind(ev.polygon, tmp)
-## Now draw the eigenvalues and the approximations
-pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA) +
-geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-## and cut the plot
-if (only.bulk == FALSE) {
-pl <- pl + coord_cartesian(xlim = 1.25 * range(Re(Approx$M.eigenvalues)),
-ylim = 1.25 * range(Im(Approx$M.eigenvalues)))
-}
-pl
-only.bulk <- TRUE
-## Check if the eigenvalues of M have been computed already
-if (is.null(Approx$M.eigenvalues)){
-Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
-}
-## Data frame for the eigenvalues of M
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-## Ellipse and extra eigenvalue for Tang et al. criteria
-st <- Approx$TangEtAl.stats
-radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
-radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius.h, radius.v)
-tmp$Type <- "Tang et al."
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "Tang et al."))
-## For cutting the graph
-bulkRe <- range(ev.polygon$Real)
-bulkIm <- range(ev.polygon$Imaginary)
-## Circle and Ellipse for eyeball approximation
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(st$center.A, st$radius.A, st$radius.A)
-tmp$Type <- "eyeball, A"
-ev.polygon <- rbind(ev.polygon, tmp)
-tmp <- eye.ellipse.df(st$radius.A + st$center.A, st$radius.B.h, st$radius.B.v)
-tmp$Type <- "eyeball, shifted B"
-ev.polygon <- rbind(ev.polygon, tmp)
-## Now draw the eigenvalues and the approximations
-pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA) +
-geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-## and cut the plot
-if (only.bulk == FALSE) {
-pl <- pl + coord_cartesian(xlim = 1.25 * range(Re(Approx$M.eigenvalues)),
-ylim = 1.25 * range(Im(Approx$M.eigenvalues)))
-}
-else {
-## Cut according to the polygons
-pl <- pl + coord_cartesian(xlim = 1.05 * range(ev.polygon$bulkRe),
-ylim = 1.05 * range(ev.polygon$bulkIm))
-}
-## Check if the eigenvalues of M have been computed already
-if (is.null(Approx$M.eigenvalues)){
-Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
-}
-## Data frame for the eigenvalues of M
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-## Ellipse and extra eigenvalue for Tang et al. criteria
-st <- Approx$TangEtAl.stats
-radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
-radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius.h, radius.v)
-tmp$Type <- "Tang et al."
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "Tang et al."))
-## For cutting the graph
-bulkRe <- range(ev.polygon$Real)
-bulkIm <- range(ev.polygon$Imaginary)
-## Circle and Ellipse for eyeball approximation
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(st$center.A, st$radius.A, st$radius.A)
-tmp$Type <- "eyeball, A"
-ev.polygon <- rbind(ev.polygon, tmp)
-tmp <- eye.ellipse.df(st$radius.A + st$center.A, st$radius.B.h, st$radius.B.v)
-tmp$Type <- "eyeball, shifted B"
-ev.polygon <- rbind(ev.polygon, tmp)
-## Now draw the eigenvalues and the approximations
-pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA) +
-geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-## and cut the plot
-if (only.bulk == FALSE) {
-pl <- pl + coord_cartesian(xlim = 1.25 * range(Re(Approx$M.eigenvalues)),
-ylim = 1.25 * range(Im(Approx$M.eigenvalues)))
-} else {
-## Cut according to the polygons
-pl <- pl + coord_cartesian(xlim = 1.05 * range(ev.polygon$bulkRe),
-ylim = 1.05 * range(ev.polygon$bulkIm))
-}
-bulkIm
-bulkRe
-## Check if the eigenvalues of M have been computed already
-if (is.null(Approx$M.eigenvalues)){
-Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
-}
-## Data frame for the eigenvalues of M
-ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
-pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
-## Data frame for approximations
-ev.polygon <- data.frame()
-ev.extra <- data.frame()
-## Circle and extra eigenvalue for May's stability criterion
-st <- Approx$May.stats
-radius <- sqrt((st$S - 1) * st$V)
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius, radius)
-tmp$Type <- "May"
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
-## Ellipse and extra eigenvalue for Tang et al. criteria
-st <- Approx$TangEtAl.stats
-radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
-radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
-center <- st$d - st$E
-tmp <- eye.ellipse.df(center, radius.h, radius.v)
-tmp$Type <- "Tang et al."
-ev.polygon <- rbind(ev.polygon, tmp)
-ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "Tang et al."))
-## For cutting the graph
-bulkRe <- range(ev.polygon$Real)
-bulkIm <- range(ev.polygon$Imaginary)
-## Circle and Ellipse for eyeball approximation
-st <- Approx$eyeball.stats
-tmp <- eye.ellipse.df(st$center.A, st$radius.A, st$radius.A)
-tmp$Type <- "eyeball, A"
-ev.polygon <- rbind(ev.polygon, tmp)
-tmp <- eye.ellipse.df(st$radius.A + st$center.A, st$radius.B.h, st$radius.B.v)
-tmp$Type <- "eyeball, shifted B"
-ev.polygon <- rbind(ev.polygon, tmp)
-## Now draw the eigenvalues and the approximations
-pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA) +
-geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 1.25)
-## and cut the plot
-if (only.bulk == FALSE) {
-pl <- pl + coord_cartesian(xlim = 1.25 * range(Re(Approx$M.eigenvalues)),
-ylim = 1.25 * range(Im(Approx$M.eigenvalues)))
-} else {
-## Cut according to the polygons
-pl <- pl + coord_cartesian(xlim = 1.05 * bulkRe,
-ylim = 1.05 * bulkIm)
-}
-pl
-devtools::load_all(".")
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner")))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade")))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade")), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.1)), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.1, mux = 0.1)), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.1, mux = -0.1)), only.bulk = TRUE)
-devtools::load_all(".")
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.1, mux = -0.1)), only.bulk = TRUE)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.1, mux = -0.1)), only.bulk = TRUE)
-devtools::load_all(".")
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.09, mux = -0.1)), only.bulk = TRUE)
-library(eyeball)
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.09, mux = -0.1)), only.bulk = TRUE)
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.09, mux = -0.1)), only.bulk = TRUE)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.09, mux = -0.1)), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1000, distribution.pairs = "FourCorner", foodweb.model = "Cascade", muy = 0.07, mux = -0.1)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 1)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 1, foodweb.model = "Cascade")))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 1, foodweb.model = "Cascade")), only.bulk = TRUE)
-devtools::load_all(".")
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 1, foodweb.model = "Cascade")), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 1, foodweb.model = "Cascade")), only.bulk = FALSE)
-devtools::load_all(".")
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 1)), only.bulk = FALSE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.25)), only.bulk = FALSE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = FALSE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-devtools::load_all(".")
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-devtools::load_all(".")
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15)), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, muy = 5)), only.bulk = TRUE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, muy = 5)), only.bulk = FALSE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 23)), only.bulk = FALSE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 23, sigmay = 48)), only.bulk = FALSE)
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 23, sigmay = 48)), only.bulk = FALSE)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 2, sigmay = 1)), only.bulk = FALSE)
-eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 2, sigmay = 1))
-eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 2, sigmay = 1)
-M <- eye.buildfoodweb.and.parameterize.M(S = 500, C = 0.15, sigmax = 2, sigmay = 1)
-eye.approximate.ReL1(M)
-eye.plot.approximations(eye.approximate.ReL1(M))
-eye.plot.approximations(eye.approximate.ReL1(M), only.bulk = TRUE)
-?eye.foodweb.cascade
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M()))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -5)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)))
-M  <- eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)
-diag(M)
-eye.approximate.ReL1(M)
-devtools::load_all(".")
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)))
-AA <- (eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)))
-AA$ReL1.observed
-AA$ReL1.May
-AA$ReL1.TangEtAl
-AA$ReL1.eyeball
-library(eyeball)
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)))
-library(eyeball)
-library(eyeball)
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)))
-devtools::load_all(".")
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50, sigma.diagonal = 2)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50, sigma.diagonal = 1)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50, sigma.diagonal = 5)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50, sigma.diagonal = 10)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50, sigma.diagonal = 10)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 600, mu.diagonal = -50, sigma.diagonal = 0.1)))
-eye.plot.approximations(eye.approximate.ReL1(eye.buildfoodweb.and.parameterize.M(S = 1600, mu.diagonal = -50, sigma.diagonal = 0.1)))
-library(eyeball)
-devtools::load_all(".")
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model =
-)
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model =
-)
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade")
-eye.plot.approximations(M1)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1000)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1000, C = 1)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1000, C = 1, distribution.pairs = "Normal", mux = -1, muy = 0.5)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1000, C = 1, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 0.5, sigmay = 0.5)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1000, C = 0.5, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 0.5, sigmay = 0.5)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1000, C = 0.5, distribution.pairs = "Normal", mux = -1, muy = 0.5, sigmax = 1, sigmay = 1)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1500, C = 0.5, distribution.pairs = "Normal", mux = -1, muy = 0.75, sigmax = 1, sigmay = 1)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1500, C = 0.5, distribution.pairs = "Normal", mux = -1, muy = 0.75, sigmax = 1, sigmay = 1, rhoxy = 0)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1500, C = 0.5, distribution.pairs = "Normal", mux = -1, muy = 0.75, sigmax = 0.5, sigmay = 0.5, rhoxy = 0)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-M1 <- eye.buildfoodweb.and.parameterize.M(foodweb.model = "Cascade", S = 1500, C = 0.5, distribution.pairs = "Normal", mux = -3, muy = 2.75, sigmax = 2.5, sigmay = 2.5, rhoxy = 0)
-eye.plot.approximations(eye.approximate.ReL1(M1))
-library(eyeball)
-devtools::load_all(".")
-devtools::load_all(".")
-?eye.foodweb.cascade
-devtools::load_all(".")
-eye.foodweb.cascade()
-devtools::load_all(".")
-library(eyeball)
-library(eyeball)
-library(eyeball)
-library(eyeball)
-library(eyeball)
-library(eyeball)
-library(eyeball)
-library(eyeball)
diff --git a/.Rproj.user/5B9FC3CC/pcs/source-pane.pper b/.Rproj.user/5B9FC3CC/pcs/source-pane.pper
index d3d70fa..1743e40 100644
--- a/.Rproj.user/5B9FC3CC/pcs/source-pane.pper
+++ b/.Rproj.user/5B9FC3CC/pcs/source-pane.pper
@@ -1,3 +1,3 @@
 {
-    "activeTab" : 6
+    "activeTab" : 0
 }
\ No newline at end of file
diff --git a/.Rproj.user/5B9FC3CC/pcs/windowlayoutstate.pper b/.Rproj.user/5B9FC3CC/pcs/windowlayoutstate.pper
index f91527d..e870d27 100644
--- a/.Rproj.user/5B9FC3CC/pcs/windowlayoutstate.pper
+++ b/.Rproj.user/5B9FC3CC/pcs/windowlayoutstate.pper
@@ -1,14 +1,14 @@
 {
     "left" : {
-        "panelheight" : 1085,
-        "splitterpos" : 629,
+        "panelheight" : 1113,
+        "splitterpos" : 841,
         "topwindowstate" : "NORMAL",
-        "windowheight" : 1123
+        "windowheight" : 1151
     },
     "right" : {
-        "panelheight" : 1085,
-        "splitterpos" : 464,
+        "panelheight" : 1113,
+        "splitterpos" : 475,
         "topwindowstate" : "NORMAL",
-        "windowheight" : 1123
+        "windowheight" : 1151
     }
 }
\ No newline at end of file
diff --git a/.Rproj.user/5B9FC3CC/pcs/workbench-pane.pper b/.Rproj.user/5B9FC3CC/pcs/workbench-pane.pper
index 1e279ed..bc9a124 100644
--- a/.Rproj.user/5B9FC3CC/pcs/workbench-pane.pper
+++ b/.Rproj.user/5B9FC3CC/pcs/workbench-pane.pper
@@ -1,4 +1,4 @@
 {
     "TabSet1" : 2,
-    "TabSet2" : 0
+    "TabSet2" : 1
 }
\ No newline at end of file
diff --git a/.Rproj.user/5B9FC3CC/persistent-state b/.Rproj.user/5B9FC3CC/persistent-state
index 5064d3c..912517f 100644
--- a/.Rproj.user/5B9FC3CC/persistent-state
+++ b/.Rproj.user/5B9FC3CC/persistent-state
@@ -1,6 +1,6 @@
 build-last-errors="[]"
 build-last-errors-base-dir="~/RPackages/eyeball/"
-build-last-outputs="[{\"output\":\"==> R CMD INSTALL --no-multiarch --with-keep.source eyeball\\n\\n\",\"type\":0},{\"output\":\"* installing to library ‘/home/sallesina/R/x86_64-pc-linux-gnu-library/3.1’\\n\",\"type\":1},{\"output\":\"* installing *source* package ‘eyeball’ ...\\n\",\"type\":1},{\"output\":\"** R\\n\",\"type\":1},{\"output\":\"** preparing package for lazy loading\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** help\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"*** installing help indices\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** building package indices\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** testing if installed package can be loaded\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"* DONE (eyeball)\\n\",\"type\":1},{\"output\":\"\",\"type\":1}]"
+build-last-outputs="[{\"output\":\"==> R CMD INSTALL --no-multiarch --with-keep.source eyeball\\n\\n\",\"type\":0},{\"output\":\"* installing to library ‘/home/sallesina/R/x86_64-pc-linux-gnu-library/3.1’\\n\",\"type\":1},{\"output\":\"* installing *source* package ‘eyeball’ ...\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** R\\n\",\"type\":1},{\"output\":\"** preparing package for lazy loading\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** help\\n\",\"type\":1},{\"output\":\"*** installing help indices\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** building package indices\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"** testing if installed package can be loaded\\n\",\"type\":1},{\"output\":\"\",\"type\":1},{\"output\":\"* DONE (eyeball)\\n\",\"type\":1},{\"output\":\"\",\"type\":1}]"
 compile_pdf_state="{\"errors\":[],\"output\":\"\",\"running\":false,\"tab_visible\":false,\"target_file\":\"\"}"
 console_procs="[]"
 files.monitored-path=""
diff --git a/.gitignore b/.gitignore
index f390c9d..c62ed59 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,5 +1,4 @@
 *~
-Code/FigureComparison/Matrices
-Code/SortedComparison/AdjacencyMatrices
-
 .Rproj.user
+.Rhistory
+
diff --git a/DESCRIPTION b/DESCRIPTION
index 7942000..c2bc637 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: eyeball
 Version: 0.1
-Date: 2015-06-15
+Date: 2015-06-19
 Title: Approximations for the real part of the leading eigenvalue of large
     community matrices.
 Description: This package provides a set of functions to: a) Build food webs b)
diff --git a/NAMESPACE b/NAMESPACE
index 71ce8b0..075e6d6 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -6,6 +6,7 @@ export(eye.ellipse.df)
 export(eye.foodweb.cascade)
 export(eye.foodweb.file)
 export(eye.foodweb.niche)
+export(eye.foodweb.random)
 export(eye.pairs.from.empirical)
 export(eye.pairs.from.fourcorner)
 export(eye.pairs.from.normal)
diff --git a/R/.Rhistory b/R/.Rhistory
deleted file mode 100644
index e69de29..0000000
diff --git a/R/eye.approximations.R b/R/eye.approximations.R
index 7754006..1266064 100644
--- a/R/eye.approximations.R
+++ b/R/eye.approximations.R
@@ -1,6 +1,35 @@
+#' Calculates various approximations for the real part of the leading eigenvalue
+#' of the community matrix M.
+#'
+#' It returns the expectation for the real part of the rightmost eigenvalue of M
+#' according to May (Nature 1972), Tang et al (Ecology Letters 2014), and
+#' Allesina et al. (Nature Communications 2015) derivations.
+#'
+#'
+#' @param M The community matrix. Note that for the Allesina et al. (eyeball) approximation, M should contain (almost) no
+#'        positive coefficients in the upper triangular part and (almost) no negative coefficients in the lower-triangular
+#'        part
+#' @param calculate_eigenvalues (TRUE by default) calculate the actual eigenvalues of M and return them
+#'
+#' @return A list composed of:
+#' \describe{
+#' \item{$M}{The community matrix}
+#' \item{$M.eigenvalues}{The eigenvalues of M (if calculate_eigenvalues = TRUE)}
+#' \item{$ReL1.observed }{The actual real part of the rightmost eigenvalue of M}
+#' \item{$ReL1.May }{The prediction according to May's approximation}
+#' \item{$ReL1.TangEtAl}{The prediction according to Tang et al. approximation}
+#' \item{$ReL1.eyeball}{The prediction according to the "eyeball" approach in Allesina et al. }
+#' \item{$May.stats}{List containing the relevant statistics for computing May's approximation}
+#' \item{$TangEtAl.stats}{List containing the relevant statistics for computing the approximation of Tang et al.}
+#' \item{$eyeball.stats}{List containing the relevant statistics for computing the eyeball approximation of Allesina et al.}
+#' }
+#'
+#' @examples
+#' M <- eye.parameterize.M(eye.foodweb.cascade())
+#' eye.approximate.ReL1(M)
 #' @export
-eye.approximate.ReL1 <- function(M, calculate.eigenvalues = TRUE){
-  if (calculate.eigenvalues == TRUE){
+eye.approximate.ReL1 <- function(M, calculate_eigenvalues = TRUE){
+  if (calculate_eigenvalues == TRUE){
     ev <- eigen(M, only.values = TRUE, symmetric = FALSE)$values
   }
   S <- dim(M)[1]
@@ -46,16 +75,9 @@ eye.approximate.ReL1 <- function(M, calculate.eigenvalues = TRUE){
   ReL1.observed <- NULL
   M.eigenvalues <- NULL
   ## if desired, calculate the actual ReL1.observed
-  if (calculate.eigenvalues == TRUE){
+  if (calculate_eigenvalues == TRUE){
     ReL1.observed <- max(Re(ev))
     M.eigenvalues <- ev
-    ## Debug
-    #print(plot(ev, xlim = c(min(Re(ev)), 1 + max(c(ReL1.observed, ReL1.May, ReL1.TangEtAl)))))
-    #print(abline(v = ReL1.observed, col = "black"))
-    #print(abline(v = ReL1.May, col = "red"))
-    #print(abline(v = ReL1.TangEtAl, col = "blue"))
-    #print(abline(v = ReL1.eyeball, col = "pink"))
-    ## End Debug
   }
   return(list(M = M,
               M.eigenvalues = M.eigenvalues,
diff --git a/R/eye.buildcommunitymatrix.R b/R/eye.buildcommunitymatrix.R
index cb7f7d7..c140176 100644
--- a/R/eye.buildcommunitymatrix.R
+++ b/R/eye.buildcommunitymatrix.R
@@ -1,32 +1,48 @@
-#' @export
+#' Builds a community matrix M based on a food web, and a parameterization for the bivariate distribution
+#'
+#' @param FW A food web object, produced by the method eye.foodweb.cascade or similar
+#' @param distribution_pairs A bivariate distribution for the pairs (possibilities are "Normal", "FourCorner", or an empirical distribution, represented by a Nx2 matrix of pairs)
+#' @param mux Mean of the negative interaction strengths
+#' @param muy Mean of the positive interaction strengths
+#' @param sigmax Standard deviation of the negative interaction strengths
+#' @param sigmay Standard deviation of the positive interaction strengths
+#' @param rhoxy Correlation between the interaction strengths
+#' @param mu_diagonal (0 by default) mean of the diagonal elements
+#' @param sigma_diagonal (0 by default) standard deviation for the diagonal elements
 #'
+#' @return The community matrix M
+#' @examples
+#' M <- eye.parameterize.M(eye.foodweb.cascade())
+#' M <- eye.parameterize.M(eye.foodweb.niche(size = 30, connectance = 0.2),
+#'      distribution_pairs = "FourCorner", mux = -1, muy = 1, sigmax = 2, sigmay = 1, rhoxy = -0.5)
+#' @export
 eye.parameterize.M <- function(FW,
-                               distribution.pairs = "Normal",
+                               distribution_pairs = "Normal",
                                mux = -1,
                                muy = 0.5,
                                sigmax = 1/4,
                                sigmay = 1/4,
                                rhoxy = -2/3,
-                               mu.diagonal = 0,
-                               sigma.diagonal = 0){
-  S <- FW$S
+                               mu_diagonal = 0,
+                               sigma_diagonal = 0){
+  S <- FW$size
   M <- matrix(0, S, S)
   ## Get off-diagonal coefficients
   NumPairs <- FW$L
-  if (is.matrix(distribution.pairs) == TRUE){
-    my.pairs <- eye.pairs.from.empirical(NumPairs = NumPairs,
-                                         Empirical.Distribution = distribution.pairs)
+  if (is.matrix(distribution_pairs) == TRUE){
+    my.pairs <- eye.pairs.from.empirical(num_pairs = NumPairs,
+                                         empirical_distribution = distribution_pairs)
   } else {
-    if (distribution.pairs == "Normal"){
-      my.pairs <- eye.pairs.from.normal(NumPairs = NumPairs,
+    if (distribution_pairs == "Normal"){
+      my.pairs <- eye.pairs.from.normal(num_pairs = NumPairs,
                                         mux = mux,
                                         muy = muy,
                                         sigmax = sigmax,
                                         sigmay = sigmay,
                                         rhoxy = rhoxy)
     }
-    if (distribution.pairs == "FourCorner"){
-      my.pairs <- eye.pairs.from.fourcorner(NumPairs = NumPairs,
+    if (distribution_pairs == "FourCorner"){
+      my.pairs <- eye.pairs.from.fourcorner(num_pairs = NumPairs,
                                           mux = mux,
                                           muy = muy,
                                           sigmax = sigmax,
@@ -37,39 +53,58 @@ eye.parameterize.M <- function(FW,
   M[FW$links] <- my.pairs[ , 1]
   M[FW$links[,2:1]] <- my.pairs[ , 2]
   ## Optional: set diagonal
-  diag(M) <- rnorm(S, mean = mu.diagonal, sd = sigma.diagonal)
+  diag(M) <- rnorm(S, mean = mu_diagonal, sd = sigma_diagonal)
   return(M)
 }
 
+#' Builds a community matrix M based on a food web model, and a parameterization for the bivariate distribution
+#'
+#' @param foodweb_model A food web model. Possible options are "Random", "Cascade", "Niche" or a file name containing an adjacency matrix
+#' @param size The number of species in the foodweb
+#' @param connectance The connectance of the food web
+#' @param distribution_pairs A bivariate distribution for the pairs (possibilities are "Normal", "FourCorner", or an empirical distribution, represented by a Nx2 matrix of pairs)
+#' @param mux Mean of the negative interaction strengths
+#' @param muy Mean of the positive interaction strengths
+#' @param sigmax Standard deviation of the negative interaction strengths
+#' @param sigmay Standard deviation of the positive interaction strengths
+#' @param rhoxy Correlation between the interaction strengths
+#' @param mu_diagonal (0 by default) mean of the diagonal elements
+#' @param sigma_diagonal (0 by default) standard deviation for the diagonal elements
+#'
+#' @return The community matrix M
+#' @examples
+#' M <- eye.buildfoodweb.and.parameterize.M(foodweb_model = "Cascade", size = 100, connectance = 0.2,
+#'      distribution_pairs = "Normal", mux = -1, muy = 1, sigmax = 1, sigmay = 1, rhoxy = -0.2,
+#'      mu_diagonal = -3, sigma_diagonal = 0)
 #' @export
-#' @import igraph
-#' @import MASS
-#' @import ggplot2
-eye.buildfoodweb.and.parameterize.M <- function(foodweb.model = "Cascade",
-                                                S = 100,
-                                                C = 0.1,
-                                                distribution.pairs = "Normal",
+eye.buildfoodweb.and.parameterize.M <- function(foodweb_model = "Cascade",
+                                                size = 100,
+                                                connectance = 0.1,
+                                                distribution_pairs = "Normal",
                                                 mux = -3,
                                                 muy = 1.5,
                                                 sigmax = 1,
                                                 sigmay = 0.75,
                                                 rhoxy = -2/3,
-                                                mu.diagonal = 0,
-                                                sigma.diagonal = 0){
+                                                mu_diagonal = 0,
+                                                sigma_diagonal = 0){
   FW <- NULL
   # if it's a file
-  if (file.exists(foodweb.model) == TRUE){
-    FW <- eye.foodweb.file(foodweb.model)
+  if (file.exists(foodweb_model) == TRUE){
+    FW <- eye.foodweb.file(foodweb_model)
+  }
+  if (foodweb_model == "Cascade"){
+    FW <- eye.foodweb.cascade(size, connectance)
   }
-  if (foodweb.model == "Cascade"){
-    FW <- eye.foodweb.cascade(S, C)
+  if (foodweb_model == "Niche"){
+    FW <- eye.foodweb.niche(size, connectance)
   }
-  if (foodweb.model == "Niche"){
-    FW <- eye.foodweb.niche(S, C)
+  if (foodweb_model == "Random"){
+    FW <- eye.foodweb.niche(size, connectance)
   }
   if (is.null(FW)){
-    stop("Invalid food web model. Please enter Cascade, Niche, or a file name for the adjacency matrix.")
+    stop("Invalid food web model. Please enter Random, Cascade, Niche, or a file name for the adjacency matrix.")
   }
-  M <- eye.parameterize.M(FW, distribution.pairs, mux, muy, sigmax, sigmay, rhoxy, mu.diagonal, sigma.diagonal)
+  M <- eye.parameterize.M(FW, distribution_pairs, mux, muy, sigmax, sigmay, rhoxy, mu_diagonal, sigma_diagonal)
   return(M)
 }
diff --git a/R/eye.buildfoodwebs.R b/R/eye.buildfoodwebs.R
index 4e568b2..f892068 100644
--- a/R/eye.buildfoodwebs.R
+++ b/R/eye.buildfoodwebs.R
@@ -1,22 +1,90 @@
+#' Builds a food web according to the Erdos-Renyi random graph model
+#'
+#' Any two species are connected with probability connectance
+#' @param size The number of species (integer)
+#' @param connectance The desired connectance (between 0 and 1)
+#'
+#'
+#' @return A connected (one piece) food web of connectance close to the desired value.
+#' The food web is a list containing:
+#' \describe{
+#' \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
+#' \item{$size}{The number of species}
+#' \item{$L}{The number of links}
+#' \item{$connectance}{The connectance L / (size choose 2)}
+#' \item{$model}{A code for the model used to generate the food web}
+#' \item{$cycles}{Boolean, whether cycles are present}
+#' }
+#' All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
+#' no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
+#' of the corresponding adjacency matrix.
+#'
+#' If the resulting network has connectance too different from the desired one, or it is not connected,
+#' the function will call itself recursively. If too many recursions are performed, it will fail.
+#' @examples
+#' FW <- eye.foodweb.random() # size = 100, connectance = 0.25
+#' FW <- eye.foodweb.random(30, 0.2)
+#' FW <- eye.foodweb.random(size = 250, connectance = 0.314)
+#' @export
+#' @import igraph
+eye.foodweb.random <- function(size = 100, connectance = 0.25){
+  ## Build the adjacency matrix
+  K <- (matrix(runif(size * size), size, size) < connectance * size / (2 * (size - 1))) * 1
+  diag(K) <- 0
+  ## Now check that it passes all the following criteria
+  ## a) the number of connections is about right
+  L <- sum(K)
+  expected.L <- choose(size, 2) * connectance
+  variance.L <- choose(size, 2) * connectance * (1. - connectance)
+  if (L < (expected.L - 2 * sqrt(variance.L)) |
+      L > (expected.L + 2 * sqrt(variance.L))){
+    warning("Failed to build foodweb, wrong number of links. Trying again...")
+    return(eye.foodweb.random(size, connectance))
+  }
+  ## b) The network is connected
+  g <- graph.adjacency(K, mode = "directed")
+  if (is.connected(g, mode = "weak") == FALSE){
+    warning("Failed to build foodweb, not connected. Trying again...")
+    return(eye.foodweb.random(size, connectance))
+  }
+  ## c) Check if it has no cycles. If so, perform a topological sort
+  ts <- topological.sort(g)
+  if (length(ts) == size){
+    K <- K[ts, ts]
+    g <- graph.adjacency(K, mode = "directed")
+    Cycles <- FALSE
+  } else {
+    Cycles <- TRUE
+  }
+  ## If it passes the tests, then build the foodweb
+  fw <- list()
+  fw$links <- get.edgelist(g)
+  fw$size <- size
+  fw$L <- dim(fw$links)[1]
+  fw$connectance <- fw$L / (size * (size - 1) * 0.5)
+  fw$model <- "Random"
+  fw$cycles <- Cycles
+  return(fw)
+}
+
 #' Builds a food web according to the cascade model
 #'
 #' The cascade model is described in Cohen et al. "Community food webs: Data and theory", 1990.
-#' The species are sorted in order from 1...S and each species has probability C of consuming any of the
+#' The species are sorted in order from 1...size and each species has probability connectance of consuming any of the
 #' preceding species.
-#' @param S The number of species (integer)
-#' @param C The desired connectance (between 0 and 1)
+#' @param size The number of species (integer)
+#' @param connectance The desired connectance (between 0 and 1)
 #'
-#' @export
 #'
 #' @return A connected (one piece) food web of connectance close to the desired value.
 #' The food web is a list containing:
 #' \describe{
 #' \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
-#' \item{$S}{The number of species}
+#' \item{$size}{The number of species}
 #' \item{$L}{The number of links}
-#' \item{$C}{The connectance L / (S choose 2)}
-#' \item{$Model}{A code for the model used to generate the food web}
-#' \item{$Cycles}{Boolean, whether cycles are present}
+#' \item{$connectance}{The connectance L / (size choose 2)}
+#' \item{$model}{A code for the model used to generate the food web}
+#' \item{$cycles}{Boolean, whether cycles are present}
 #' }
 #' All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
 #' no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
@@ -25,60 +93,61 @@
 #' If the resulting network has connectance too different from the desired one, or it is not connected,
 #' the function will call itself recursively. If too many recursions are performed, it will fail.
 #' @examples
-#' FW <- eye.foodweb.cascade() # S = 100, C = 0.25
+#' FW <- eye.foodweb.cascade() # size = 100, connectance = 0.25
 #' FW <- eye.foodweb.cascade(30, 0.2)
-#' FW <- eye.foodweb.cascade(S = 250, C = 0.314)
-eye.foodweb.cascade <- function(S = 100, C = 0.25){
+#' FW <- eye.foodweb.cascade(size = 250, connectance = 0.314)
+#' @export
+#' @import igraph
+eye.foodweb.cascade <- function(size = 100, connectance = 0.25){
   ## Build the adjacency matrix
-  K <- (matrix(runif(S * S), S, S) < C) * 1
+  K <- (matrix(runif(size * size), size, size) < connectance) * 1
   K[lower.tri(K)] <- 0
   diag(K) <- 0
   ## Now check that it passes all the following criteria
   ## a) the number of connections is about right
   L <- sum(K)
-  expected.L <- choose(S, 2) * C
-  variance.L <- choose(S, 2) * C * (1. - C)
+  expected.L <- choose(size, 2) * connectance
+  variance.L <- choose(size, 2) * connectance * (1. - connectance)
   if (L < (expected.L - 2 * sqrt(variance.L)) |
       L > (expected.L + 2 * sqrt(variance.L))){
         warning("Failed to build foodweb, wrong number of links. Trying again...")
-        return(eye.foodweb.cascade(S, C))
+        return(eye.foodweb.cascade(size, connectance))
       }
   ## b) The network is connected
   g <- graph.adjacency(K, mode = "directed")
   if (is.connected(g, mode = "weak") == FALSE){
     warning("Failed to build foodweb, not connected. Trying again...")
-    return(eye.foodweb.cascade(S, C))
+    return(eye.foodweb.cascade(size, connectance))
   }
   ## If it passes the tests, then build the foodweb
   fw <- list()
   fw$links <- get.edgelist(g)
-  fw$S <- S
+  fw$size <- size
   fw$L <- dim(fw$links)[1]
-  fw$C <- fw$L / (S * (S - 1) * 0.5)
-  fw$Model <- "Cascade"
-  fw$Cycles <- FALSE
+  fw$connectance <- fw$L / (size * (size - 1) * 0.5)
+  fw$model <- "Cascade"
+  fw$cycles <- FALSE
   return(fw)
 }
 
 #' Builds a food web according to the niche model
 #'
 #' The niche model is described in Williams and Martinez, Nature 2000.
-#' The species are sorted in order from 1...S and each species consumes a set of
+#' The species are sorted in order from 1...size and each species consumes a set of
 #' adjacent species.
-#' @param S The number of species (integer)
-#' @param C The desired connectance (between 0 and 1)
+#' @param size The number of species (integer)
+#' @param connectance The desired connectance (between 0 and 1)
 #'
-#' @export
 #'
 #' @return A food web
 #' The food web is a list containing:
 #' \describe{
 #' \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
-#' \item{$S}{The number of species}
+#' \item{$size}{The number of species}
 #' \item{$L}{The number of links}
-#' \item{$C}{The connectance L / (S choose 2)}
-#' \item{$Model}{A code for the model used to generate the food web}
-#' \item{$Cycles}{Boolean, whether cycles are present}
+#' \item{$connectance}{The connectance L / (size choose 2)}
+#' \item{$model}{A code for the model used to generate the food web}
+#' \item{$cycles}{Boolean, whether cycles are present}
 #' }
 #' All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
 #' no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
@@ -87,25 +156,27 @@ eye.foodweb.cascade <- function(S = 100, C = 0.25){
 #' If the resulting network has connectance too different from the desired one, or it is not connected,
 #' the function will call itself recursively. If too many recursions are performed, it will fail.
 #' @examples
-#' FW <- eye.foodweb.niche() # S = 100, C = 0.25
+#' FW <- eye.foodweb.niche() # size = 100, connectance = 0.25
 #' FW <- eye.foodweb.niche(30, 0.2)
-#' FW <- eye.foodweb.niche(S = 250, C = 0.314)
-eye.foodweb.niche <- function(S = 100, C = 0.25){
+#' FW <- eye.foodweb.niche(size = 250, connectance = 0.314)
+#' @export
+#' @import igraph
+eye.foodweb.niche <- function(size = 100, connectance = 0.25){
   ## Build the adjacency matrix
-  K <- matrix(0, S, S)
+  K <- matrix(0, size, size)
   ## assign the niche value for all species
   ## Note: the niche values of the species are ordered
-  ni <- sort(runif(S, 0, 1))
+  ni <- sort(runif(size, 0, 1))
   ## determine the radius of the food spectrum for each species
   ## the radius is drawn randomly from a ni * Beta(1, beta) distribution
-  beta <- 1.0 / C - 1
-  ri <- rbeta(S, 1, beta) * ni
+  beta <- 1.0 / connectance - 1
+  ri <- rbeta(size, 1, beta) * ni
   ## Set species with the lowest niche value to be the basal species
   ri[1] <- 0
   ## The center of the food spectrum for each species are
   ## drawn uniformly from an uniform distribution on [ri/2, ni]
-  ci <- numeric(S)
-  for(i in 1:S){
+  ci <- numeric(size)
+  for(i in 1:size){
     ci[i] <- runif(1, ri[i] / 2, min(ni[i], 1 - ri[i] / 2))
     ## Correction as in Allesina et al Science 2008
     ## Determine the boundary for the food spectrum for each species:
@@ -114,7 +185,7 @@ eye.foodweb.niche <- function(S = 100, C = 0.25){
     ## check which species are falling into that interval: [lower, upper]
     ## and set the corresponding value in the ajacency matrix M to be 1,
     ## indicating a link is established between the two species.
-    for(j in 1:S){
+    for(j in 1:size){
       if(ni[j] > lower & ni[j] < upper)
         K[j, i] <- 1
         ## remove double edges
@@ -132,22 +203,22 @@ eye.foodweb.niche <- function(S = 100, C = 0.25){
   ## Now check that it passes all the following criteria
   ## a) the number of connections is about right
   L <- sum(K)
-  expected.L <- choose(S, 2) * C
-  variance.L <- choose(S, 2) * C * (1. - C)
+  expected.L <- choose(size, 2) * connectance
+  variance.L <- choose(size, 2) * connectance * (1. - connectance)
   if (L < (expected.L - 2 * sqrt(variance.L)) |
         L > (expected.L + 2 * sqrt(variance.L))){
     warning("Failed to build foodweb, wrong number of links. Trying again...")
-    return(eye.foodweb.niche(S, C))
+    return(eye.foodweb.niche(size, connectance))
   }
   ## b) The network is connected
   g <- graph.adjacency(K, mode = "directed")
   if (is.connected(g, mode = "weak") == FALSE){
     warning("Failed to build foodweb, not connected. Trying again...")
-    return(eye.foodweb.niche(S, C))
+    return(eye.foodweb.niche(size, connectance))
   }
   ## c) Check if it has no cycles. If so, perform a topological sort
   ts <- topological.sort(g)
-  if (length(ts) == S){
+  if (length(ts) == size){
     K <- K[ts, ts]
     g <- graph.adjacency(K, mode = "directed")
     Cycles <- FALSE
@@ -157,11 +228,11 @@ eye.foodweb.niche <- function(S = 100, C = 0.25){
   ## If it passes the tests, then build the foodweb
   fw <- list()
   fw$links <- get.edgelist(g)
-  fw$S <- S
+  fw$size <- size
   fw$L <- dim(fw$links)[1]
-  fw$C <- fw$L / (S * (S - 1) * 0.5)
-  fw$Model <- "Niche"
-  fw$Cycles <- Cycles
+  fw$connectance <- fw$L / (size * (size - 1) * 0.5)
+  fw$model <- "Niche"
+  fw$cycles <- Cycles
   return(fw)
 }
 
@@ -172,7 +243,6 @@ eye.foodweb.niche <- function(S = 100, C = 0.25){
 #' @param filename A text file containing the adjacency matrix (binary, rows = resources, cols = consumers) of a food web
 #' The food web is a list containing:
 #'
-#' @export
 #'
 #' @return
 #'
@@ -180,21 +250,23 @@ eye.foodweb.niche <- function(S = 100, C = 0.25){
 #'
 #' \describe{
 #' \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
-#' \item{$S}{The number of species}
+#' \item{$size}{The number of species}
 #' \item{$L}{The number of links}
-#' \item{$C}{The connectance L / (S choose 2)}
-#' \item{$Model}{A code for the model used to generate the food web}
-#' \item{$Cycles}{Boolean, whether cycles are present}
+#' \item{$connectance}{The connectance L / (size choose 2)}
+#' \item{$model}{A code for the model used to generate the food web}
+#' \item{$cycles}{Boolean, whether cycles are present}
 #' }
 #' All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
 #' no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
 #' of the corresponding adjacency matrix.
+#' @export
+#' @import igraph
 eye.foodweb.file <- function(filename){
   ## Build the adjacency matrix
   K <- (as.matrix(read.table(filename, header = FALSE)) > 0) * 1
   colnames(K) <- NULL
   rownames(K) <- NULL
-  S <- dim(K)[1]
+  size <- dim(K)[1]
   diag(K) <- 0
   ## Check double arrows
   KK <- K + t(K)
@@ -204,7 +276,7 @@ eye.foodweb.file <- function(filename){
   g <- graph.adjacency(K, mode = "directed")
   ## Check if it has no cycles. If so, perform a topological sort
   ts <- topological.sort(g)
-  if (length(ts) == S){
+  if (length(ts) == size){
     K <- K[ts, ts]
     g <- graph.adjacency(K, mode = "directed")
     Cycles <- FALSE
@@ -214,10 +286,10 @@ eye.foodweb.file <- function(filename){
   ## Build the foodweb
   fw <- list()
   fw$links <- get.edgelist(g)
-  fw$S <- S
+  fw$size <- size
   fw$L <- dim(fw$links)[1]
-  fw$C <- fw$L / (S * (S - 1) * 0.5)
-  fw$Model <- "From File"
-  fw$Cycles <- Cycles
+  fw$connectance <- fw$L / (size * (size - 1) * 0.5)
+  fw$model <- "From File"
+  fw$cycles <- Cycles
   return(fw)
 }
diff --git a/R/eye.plotting.R b/R/eye.plotting.R
index 3dca1c9..26ddb31 100644
--- a/R/eye.plotting.R
+++ b/R/eye.plotting.R
@@ -1,3 +1,4 @@
+#' Function needed to calculate an ellipse in the complex plane
 #' @export
 eye.ellipse.df <- function(centerx, radiusx, radiusy){
   thetas <- seq(pi / 2.0, 0.0, length = 1000)
@@ -8,20 +9,32 @@ eye.ellipse.df <- function(centerx, radiusx, radiusy){
   return(data.frame(Real = x + centerx, Imaginary = y))
 }
 
+#' Plot a graph showing the quality of the approximations
+#'
+#' @param approximation An object returned by eye.approximate.ReL1
+#' @param only_bulk If TRUE, plot only the eigenvalues close to the center of the distribution
+#'
+#' @return A ggplot plot
+#'
+#' @examples
+#' M <- eye.parameterize.M(eye.foodweb.cascade())
+#' Appr <- eye.approximate.ReL1(M)
+#' eye.plot.approximations(Appr)
 #' @export
-eye.plot.approximations <- function(Approx, only.bulk = FALSE){
+#' @import ggplot2
+eye.plot.approximations <- function(approximation, only_bulk = FALSE){
   ## Check if the eigenvalues of M have been computed already
-  if (is.null(Approx$M.eigenvalues)){
-    Approx$M.eigenvalues <- eigen(Approx$M, only.values = TRUE, symmetric = FALSE)$values
+  if (is.null(approximation$M.eigenvalues)){
+    approximation$M.eigenvalues <- eigen(approximation$M, only.values = TRUE, symmetric = FALSE)$values
   }
   ## Data frame for the eigenvalues of M
-  ev <- data.frame(Real = Re(Approx$M.eigenvalues), Imaginary = Im(Approx$M.eigenvalues))
+  ev <- data.frame(Real = Re(approximation$M.eigenvalues), Imaginary = Im(approximation$M.eigenvalues))
   pl <- ggplot(ev, aes(Real, Imaginary)) + geom_point(alpha = 0.7) + theme_bw()
   ## Data frame for approximations
   ev.polygon <- data.frame()
   ev.extra <- data.frame()
   ## Circle and extra eigenvalue for May's stability criterion
-  st <- Approx$May.stats
+  st <- approximation$May.stats
   radius <- sqrt((st$S - 1) * st$V)
   center <- st$d - st$E
   tmp <- eye.ellipse.df(center, radius, radius)
@@ -29,7 +42,7 @@ eye.plot.approximations <- function(Approx, only.bulk = FALSE){
   ev.polygon <- rbind(ev.polygon, tmp)
   ev.extra <- rbind(ev.extra, data.frame(Real = (st$S-1) * st$E + st$d, Imaginary = 0, Type = "May"))
   ## Ellipse and extra eigenvalue for Tang et al. criteria
-  st <- Approx$TangEtAl.stats
+  st <- approximation$TangEtAl.stats
   radius.h <- sqrt((st$S - 1) * st$V * (1 + st$rho))
   radius.v <- sqrt((st$S - 1) * st$V * (1 - st$rho))
   center <- st$d - st$E
@@ -41,7 +54,7 @@ eye.plot.approximations <- function(Approx, only.bulk = FALSE){
   bulkRe <- range(ev.polygon$Real)
   bulkIm <- range(ev.polygon$Imaginary)
   ## Circle and Ellipse for eyeball approximation
-  st <- Approx$eyeball.stats
+  st <- approximation$eyeball.stats
   tmp <- eye.ellipse.df(st$center.A + st$d, st$radius.A, st$radius.A)
   tmp$Type <- "eyeball, A"
   ev.polygon <- rbind(ev.polygon, tmp)
@@ -52,9 +65,9 @@ eye.plot.approximations <- function(Approx, only.bulk = FALSE){
   pl <- pl + geom_polygon(data = ev.polygon, aes(Real, Imaginary, colour = Type), fill = NA) +
     geom_point(data = ev.extra, aes(Real, Imaginary, colour = Type), shape = 1, size = 2, show_guide = FALSE)
   ## and cut the plot
-  if (only.bulk == FALSE) {
-    pl <- pl + coord_cartesian(xlim = 1.05 * range(c(bulkRe, Re(Approx$M.eigenvalues))),
-                              ylim = 1.05 * range(c(bulkIm, Im(Approx$M.eigenvalues))))
+  if (only_bulk == FALSE) {
+    pl <- pl + coord_cartesian(xlim = 1.05 * range(c(bulkRe, Re(approximation$M.eigenvalues))),
+                              ylim = 1.05 * range(c(bulkIm, Im(approximation$M.eigenvalues))))
   } else {
     ## Cut according to the polygons
     pl <- pl + coord_cartesian(xlim = 1.05 * bulkRe,
diff --git a/R/eye.sampledistributions.R b/R/eye.sampledistributions.R
index 5fa7e01..77ea02e 100644
--- a/R/eye.sampledistributions.R
+++ b/R/eye.sampledistributions.R
@@ -1,19 +1,23 @@
 #' Sample pairs of interactions from a bivariate normal distribution
 #'
-#' @param NumPairs The number of pairs to sample
+#' @param num_pairs The number of pairs to sample
 #' @param mux The desired mean for the first marginal distribution
 #' @param muy The desired mean for the second marginal distribution
 #' @param sigmax The desired standard deviation for the first marginal distribution
 #' @param sigmay The desired standard deviation for the second marginal distribution
-#' @param rhoxy The desired correlation
+#' @param rhoxy The desired correlation between the pairs of interaction strengths
 #'
-#' @export
 #'
-#' @return A NumPairs x 2 matrix of interaction strengths sampled from the distribution
+#' @return A num_pairs x 2 matrix of interaction strengths sampled from the distribution
 #' @examples
 #' mypairs <- eye.pairs.from.normal() # default values
-#' mypairs <- eye.pairs.from.normal(NumPairs = 120)
-eye.pairs.from.normal <- function(NumPairs = 10,
+#' mypairs <- eye.pairs.from.normal(450)
+#' mypairs <- eye.pairs.from.normal(num_pairs = 450)
+#' mypairs <- eye.pairs.from.normal(num_pairs = 450, mux = 0, muy = 0,
+#'            sigmax = 1, sigmay = 1, rhoxy = 0) # as for May's matrices
+#' @export
+#' @import MASS
+eye.pairs.from.normal <- function(num_pairs = 10,
                               mux = -1,
                               muy = 0.5,
                               sigmax = 1/4,
@@ -25,60 +29,63 @@ eye.pairs.from.normal <- function(NumPairs = 10,
                                 rhoxy * sigmax * sigmay,
                                 sigmay^2),
                               2, 2)
-  Pairs <- mvrnorm(NumPairs, mus, covariance.matrix)
+  Pairs <- mvrnorm(num_pairs, mus, covariance.matrix)
   return(Pairs)
 }
 
 #' Sample pairs of interactions from the "Four-Corner" distribution described
 #' in Allesina et al. 2014
 #'
-#' @param NumPairs The number of pairs to sample
+#' @param num_pairs The number of pairs to sample
 #' @param mux The desired mean for the first marginal distribution
 #' @param muy The desired mean for the second marginal distribution
 #' @param sigmax The desired standard deviation for the first marginal distribution
 #' @param sigmay The desired standard deviation for the second marginal distribution
 #' @param rhoxy The desired correlation
 #'
-#' @export
 #'
-#' @return A NumPairs x 2 matrix of interaction strengths sampled from the distribution
+#' @return A num_pairs x 2 matrix of interaction strengths sampled from the distribution
 #' @examples
 #' mypairs <- eye.pairs.from.fourcorner() # default values
-#' mypairs <- eye.pairs.from.fourcorner(NumPairs = 120)
-eye.pairs.from.fourcorner <- function(NumPairs = 10,
+#' mypairs <- eye.pairs.from.fourcorner(450)
+#' mypairs <- eye.pairs.from.fourcorner(num_pairs = 450)
+#' mypairs <- eye.pairs.from.fourcorner(num_pairs = 450, mux = 0, muy = 0,
+#'            sigmax = 1, sigmay = 1, rhoxy = 0) # as for May's matrices
+#' @export
+eye.pairs.from.fourcorner <- function(num_pairs = 10,
                                   mux = -1,
                                   muy = 0.5,
                                   sigmax = 1/4,
                                   sigmay = 1/4,
                                   rhoxy = -2/3){
   gamma <- (rhoxy + 1) / 2
-  PosNeg <- sign(rnorm(NumPairs))
-  SwitchSign <- sign(runif(NumPairs, -(1-gamma), gamma))
+  PosNeg <- sign(rnorm(num_pairs))
+  SwitchSign <- sign(runif(num_pairs, -(1-gamma), gamma))
   Pairs <- cbind(mux + PosNeg * sigmax, muy + PosNeg * SwitchSign * sigmay)
   return(Pairs)
 }
 
 #' Sample pairs of interactions from an empirical distribution
 #'
-#' @param NumPairs The number of pairs to sample
-#' @param Empirical.Distribution A Nx2 matrix with the effects of consumers on resources
+#' @param num_pairs The number of pairs to sample
+#' @param empirical_distribution A Nx2 matrix with the effects of consumers on resources
 #'                                in the first column and those of resources on consumers in the second.
 #'
-#' @export
 #'
-#' @return A NumPairs x 2 matrix of interaction strengths sampled from the distribution
+#' @return A num_pairs x 2 matrix of interaction strengths sampled from the distribution
 #' @examples
 #' mypairs <- eye.pairs.from.empirical(100, matrix(runif(20), 10, 2))
-eye.pairs.from.empirical <- function(NumPairs = 10,
-                                 Empirical.Distribution = NULL){
-  if (is.matrix(Empirical.Distribution) == TRUE){
-    NR <- dim(Empirical.Distribution)[2]
+#' @export
+eye.pairs.from.empirical <- function(num_pairs = 10,
+                                 empirical_distribution = NULL){
+  if (is.matrix(empirical_distribution) == TRUE){
+    NR <- dim(empirical_distribution)[2]
     if (NR == 2){
-      Pairs <- matrix(NumPairs, 2)
-      SampledPairs <- sample(1:NR, NumPairs, replace = TRUE)
-      Pairs <- Empirical.Distribution[SampledPairs,]
+      Pairs <- matrix(num_pairs, 2)
+      SampledPairs <- sample(1:NR, num_pairs, replace = TRUE)
+      Pairs <- empirical_distribution[SampledPairs,]
       return(Pairs)
     }
   }
-  stop("The Empirical.Distribution must be a N x 2 matrix!")
+  stop("The empirical_distribution must be a N x 2 matrix!")
 }
diff --git a/README.Md b/README.Md
new file mode 100644
index 0000000..1734275
--- /dev/null
+++ b/README.Md
@@ -0,0 +1,22 @@
+The code is provided as an R package.
+
+To install, you need the packages devtools, MASS, igraph, and ggplot2.
+You can install these packages in R with the command
+
+> install.packages(devtools)
+
+etc.
+
+To install the eyeball package, type the following in R
+
+> library(devtools)
+
+> install_github("StefanoAllesina/eyeball")
+
+> library(eyeball)
+
+That's it! For any question, comments, problem, bug, feel free to write me at
+sallesina@uchicago.edu
+
+
+
diff --git a/man/eye.approximate.ReL1.Rd b/man/eye.approximate.ReL1.Rd
new file mode 100644
index 0000000..e5582cc
--- /dev/null
+++ b/man/eye.approximate.ReL1.Rd
@@ -0,0 +1,39 @@
+% Generated by roxygen2 (4.0.1): do not edit by hand
+\name{eye.approximate.ReL1}
+\alias{eye.approximate.ReL1}
+\title{Calculates various approximations for the real part of the leading eigenvalue
+of the community matrix M.}
+\usage{
+eye.approximate.ReL1(M, calculate_eigenvalues = TRUE)
+}
+\arguments{
+\item{M}{The community matrix. Note that for the Allesina et al. (eyeball) approximation, M should contain (almost) no
+positive coefficients in the upper triangular part and (almost) no negative coefficients in the lower-triangular
+part}
+
+\item{calculate_eigenvalues}{(TRUE by default) calculate the actual eigenvalues of M and return them}
+}
+\value{
+A list composed of:
+\describe{
+\item{$M}{The community matrix}
+\item{$M.eigenvalues}{The eigenvalues of M (if calculate_eigenvalues = TRUE)}
+\item{$ReL1.observed }{The actual real part of the rightmost eigenvalue of M}
+\item{$ReL1.May }{The prediction according to May's approximation}
+\item{$ReL1.TangEtAl}{The prediction according to Tang et al. approximation}
+\item{$ReL1.eyeball}{The prediction according to the "eyeball" approach in Allesina et al. }
+\item{$May.stats}{List containing the relevant statistics for computing May's approximation}
+\item{$TangEtAl.stats}{List containing the relevant statistics for computing the approximation of Tang et al.}
+\item{$eyeball.stats}{List containing the relevant statistics for computing the eyeball approximation of Allesina et al.}
+}
+}
+\description{
+It returns the expectation for the real part of the rightmost eigenvalue of M
+according to May (Nature 1972), Tang et al (Ecology Letters 2014), and
+Allesina et al. (Nature Communications 2015) derivations.
+}
+\examples{
+M <- eye.parameterize.M(eye.foodweb.cascade())
+eye.approximate.ReL1(M)
+}
+
diff --git a/man/eye.buildfoodweb.and.parameterize.M.Rd b/man/eye.buildfoodweb.and.parameterize.M.Rd
new file mode 100644
index 0000000..ad2a17e
--- /dev/null
+++ b/man/eye.buildfoodweb.and.parameterize.M.Rd
@@ -0,0 +1,45 @@
+% Generated by roxygen2 (4.0.1): do not edit by hand
+\name{eye.buildfoodweb.and.parameterize.M}
+\alias{eye.buildfoodweb.and.parameterize.M}
+\title{Builds a community matrix M based on a food web model, and a parameterization for the bivariate distribution}
+\usage{
+eye.buildfoodweb.and.parameterize.M(foodweb_model = "Cascade", size = 100,
+  connectance = 0.1, distribution_pairs = "Normal", mux = -3, muy = 1.5,
+  sigmax = 1, sigmay = 0.75, rhoxy = -2/3, mu_diagonal = 0,
+  sigma_diagonal = 0)
+}
+\arguments{
+\item{foodweb_model}{A food web model. Possible options are "Random", "Cascade", "Niche" or a file name containing an adjacency matrix}
+
+\item{size}{The number of species in the foodweb}
+
+\item{connectance}{The connectance of the food web}
+
+\item{distribution_pairs}{A bivariate distribution for the pairs (possibilities are "Normal", "FourCorner", or an empirical distribution, represented by a Nx2 matrix of pairs)}
+
+\item{mux}{Mean of the negative interaction strengths}
+
+\item{muy}{Mean of the positive interaction strengths}
+
+\item{sigmax}{Standard deviation of the negative interaction strengths}
+
+\item{sigmay}{Standard deviation of the positive interaction strengths}
+
+\item{rhoxy}{Correlation between the interaction strengths}
+
+\item{mu_diagonal}{(0 by default) mean of the diagonal elements}
+
+\item{sigma_diagonal}{(0 by default) standard deviation for the diagonal elements}
+}
+\value{
+The community matrix M
+}
+\description{
+Builds a community matrix M based on a food web model, and a parameterization for the bivariate distribution
+}
+\examples{
+M <- eye.buildfoodweb.and.parameterize.M(foodweb_model = "Cascade", size = 100, connectance = 0.2,
+     distribution_pairs = "Normal", mux = -1, muy = 1, sigmax = 1, sigmay = 1, rhoxy = -0.2,
+     mu_diagonal = -3, sigma_diagonal = 0)
+}
+
diff --git a/man/eye.ellipse.df.Rd b/man/eye.ellipse.df.Rd
new file mode 100644
index 0000000..94f31da
--- /dev/null
+++ b/man/eye.ellipse.df.Rd
@@ -0,0 +1,11 @@
+% Generated by roxygen2 (4.0.1): do not edit by hand
+\name{eye.ellipse.df}
+\alias{eye.ellipse.df}
+\title{Function needed to calculate an ellipse in the complex plane}
+\usage{
+eye.ellipse.df(centerx, radiusx, radiusy)
+}
+\description{
+Function needed to calculate an ellipse in the complex plane
+}
+
diff --git a/man/eye.foodweb.cascade.Rd b/man/eye.foodweb.cascade.Rd
index 595ffd3..a72bf1c 100644
--- a/man/eye.foodweb.cascade.Rd
+++ b/man/eye.foodweb.cascade.Rd
@@ -3,23 +3,23 @@
 \alias{eye.foodweb.cascade}
 \title{Builds a food web according to the cascade model}
 \usage{
-eye.foodweb.cascade(S = 100, C = 0.25)
+eye.foodweb.cascade(size = 100, connectance = 0.25)
 }
 \arguments{
-\item{S}{The number of species (integer)}
+\item{size}{The number of species (integer)}
 
-\item{C}{The desired connectance (between 0 and 1)}
+\item{connectance}{The desired connectance (between 0 and 1)}
 }
 \value{
 A connected (one piece) food web of connectance close to the desired value.
 The food web is a list containing:
 \describe{
 \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
-\item{$S}{The number of species}
+\item{$size}{The number of species}
 \item{$L}{The number of links}
-\item{$C}{The connectance L / (S choose 2)}
-\item{$Model}{A code for the model used to generate the food web}
-\item{$Cycles}{Boolean, whether cycles are present}
+\item{$connectance}{The connectance L / (size choose 2)}
+\item{$model}{A code for the model used to generate the food web}
+\item{$cycles}{Boolean, whether cycles are present}
 }
 All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
 no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
@@ -30,12 +30,12 @@ the function will call itself recursively. If too many recursions are performed,
 }
 \description{
 The cascade model is described in Cohen et al. "Community food webs: Data and theory", 1990.
-The species are sorted in order from 1...S and each species has probability C of consuming any of the
+The species are sorted in order from 1...size and each species has probability connectance of consuming any of the
 preceding species.
 }
 \examples{
-FW <- eye.foodweb.cascade() # S = 100, C = 0.25
+FW <- eye.foodweb.cascade() # size = 100, connectance = 0.25
 FW <- eye.foodweb.cascade(30, 0.2)
-FW <- eye.foodweb.cascade(S = 250, C = 0.314)
+FW <- eye.foodweb.cascade(size = 250, connectance = 0.314)
 }
 
diff --git a/man/eye.foodweb.file.Rd b/man/eye.foodweb.file.Rd
index c5404aa..97538ff 100644
--- a/man/eye.foodweb.file.Rd
+++ b/man/eye.foodweb.file.Rd
@@ -14,11 +14,11 @@ The food web is a list containing:
 
 \describe{
 \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
-\item{$S}{The number of species}
+\item{$size}{The number of species}
 \item{$L}{The number of links}
-\item{$C}{The connectance L / (S choose 2)}
-\item{$Model}{A code for the model used to generate the food web}
-\item{$Cycles}{Boolean, whether cycles are present}
+\item{$connectance}{The connectance L / (size choose 2)}
+\item{$model}{A code for the model used to generate the food web}
+\item{$cycles}{Boolean, whether cycles are present}
 }
 All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
 no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
diff --git a/man/eye.foodweb.niche.Rd b/man/eye.foodweb.niche.Rd
index 16699ae..c241e0d 100644
--- a/man/eye.foodweb.niche.Rd
+++ b/man/eye.foodweb.niche.Rd
@@ -3,23 +3,23 @@
 \alias{eye.foodweb.niche}
 \title{Builds a food web according to the niche model}
 \usage{
-eye.foodweb.niche(S = 100, C = 0.25)
+eye.foodweb.niche(size = 100, connectance = 0.25)
 }
 \arguments{
-\item{S}{The number of species (integer)}
+\item{size}{The number of species (integer)}
 
-\item{C}{The desired connectance (between 0 and 1)}
+\item{connectance}{The desired connectance (between 0 and 1)}
 }
 \value{
 A food web
 The food web is a list containing:
 \describe{
 \item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
-\item{$S}{The number of species}
+\item{$size}{The number of species}
 \item{$L}{The number of links}
-\item{$C}{The connectance L / (S choose 2)}
-\item{$Model}{A code for the model used to generate the food web}
-\item{$Cycles}{Boolean, whether cycles are present}
+\item{$connectance}{The connectance L / (size choose 2)}
+\item{$model}{A code for the model used to generate the food web}
+\item{$cycles}{Boolean, whether cycles are present}
 }
 All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
 no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
@@ -30,12 +30,12 @@ the function will call itself recursively. If too many recursions are performed,
 }
 \description{
 The niche model is described in Williams and Martinez, Nature 2000.
-The species are sorted in order from 1...S and each species consumes a set of
+The species are sorted in order from 1...size and each species consumes a set of
 adjacent species.
 }
 \examples{
-FW <- eye.foodweb.niche() # S = 100, C = 0.25
+FW <- eye.foodweb.niche() # size = 100, connectance = 0.25
 FW <- eye.foodweb.niche(30, 0.2)
-FW <- eye.foodweb.niche(S = 250, C = 0.314)
+FW <- eye.foodweb.niche(size = 250, connectance = 0.314)
 }
 
diff --git a/man/eye.foodweb.random.Rd b/man/eye.foodweb.random.Rd
new file mode 100644
index 0000000..0d61d64
--- /dev/null
+++ b/man/eye.foodweb.random.Rd
@@ -0,0 +1,39 @@
+% Generated by roxygen2 (4.0.1): do not edit by hand
+\name{eye.foodweb.random}
+\alias{eye.foodweb.random}
+\title{Builds a food web according to the Erdos-Renyi random graph model}
+\usage{
+eye.foodweb.random(size = 100, connectance = 0.25)
+}
+\arguments{
+\item{size}{The number of species (integer)}
+
+\item{connectance}{The desired connectance (between 0 and 1)}
+}
+\value{
+A connected (one piece) food web of connectance close to the desired value.
+The food web is a list containing:
+\describe{
+\item{$links}{The trophic links from resource to consumer (L x 2 matrix)}
+\item{$size}{The number of species}
+\item{$L}{The number of links}
+\item{$connectance}{The connectance L / (size choose 2)}
+\item{$model}{A code for the model used to generate the food web}
+\item{$cycles}{Boolean, whether cycles are present}
+}
+All the self-loops and double arrows (a eats b, b eats a) are removed. If the food web contains
+no cycles, then it is sorted that all the coefficients would be in the upper-triangular part
+of the corresponding adjacency matrix.
+
+If the resulting network has connectance too different from the desired one, or it is not connected,
+the function will call itself recursively. If too many recursions are performed, it will fail.
+}
+\description{
+Any two species are connected with probability connectance
+}
+\examples{
+FW <- eye.foodweb.random() # size = 100, connectance = 0.25
+FW <- eye.foodweb.random(30, 0.2)
+FW <- eye.foodweb.random(size = 250, connectance = 0.314)
+}
+
diff --git a/man/eye.pairs.from.empirical.Rd b/man/eye.pairs.from.empirical.Rd
index 5aa13ae..88961a7 100644
--- a/man/eye.pairs.from.empirical.Rd
+++ b/man/eye.pairs.from.empirical.Rd
@@ -3,16 +3,16 @@
 \alias{eye.pairs.from.empirical}
 \title{Sample pairs of interactions from an empirical distribution}
 \usage{
-eye.pairs.from.empirical(NumPairs = 10, Empirical.Distribution = NULL)
+eye.pairs.from.empirical(num_pairs = 10, empirical_distribution = NULL)
 }
 \arguments{
-\item{NumPairs}{The number of pairs to sample}
+\item{num_pairs}{The number of pairs to sample}
 
-\item{Empirical.Distribution}{A Nx2 matrix with the effects of consumers on resources
+\item{empirical_distribution}{A Nx2 matrix with the effects of consumers on resources
                                in the first column and those of resources on consumers in the second.}
 }
 \value{
-A NumPairs x 2 matrix of interaction strengths sampled from the distribution
+A num_pairs x 2 matrix of interaction strengths sampled from the distribution
 }
 \description{
 Sample pairs of interactions from an empirical distribution
diff --git a/man/eye.pairs.from.fourcorner.Rd b/man/eye.pairs.from.fourcorner.Rd
index 16718ff..118569e 100644
--- a/man/eye.pairs.from.fourcorner.Rd
+++ b/man/eye.pairs.from.fourcorner.Rd
@@ -4,11 +4,11 @@
 \title{Sample pairs of interactions from the "Four-Corner" distribution described
 in Allesina et al. 2014}
 \usage{
-eye.pairs.from.fourcorner(NumPairs = 10, mux = -1, muy = 0.5,
+eye.pairs.from.fourcorner(num_pairs = 10, mux = -1, muy = 0.5,
   sigmax = 1/4, sigmay = 1/4, rhoxy = -2/3)
 }
 \arguments{
-\item{NumPairs}{The number of pairs to sample}
+\item{num_pairs}{The number of pairs to sample}
 
 \item{mux}{The desired mean for the first marginal distribution}
 
@@ -21,7 +21,7 @@ eye.pairs.from.fourcorner(NumPairs = 10, mux = -1, muy = 0.5,
 \item{rhoxy}{The desired correlation}
 }
 \value{
-A NumPairs x 2 matrix of interaction strengths sampled from the distribution
+A num_pairs x 2 matrix of interaction strengths sampled from the distribution
 }
 \description{
 Sample pairs of interactions from the "Four-Corner" distribution described
@@ -29,6 +29,9 @@ in Allesina et al. 2014
 }
 \examples{
 mypairs <- eye.pairs.from.fourcorner() # default values
-mypairs <- eye.pairs.from.fourcorner(NumPairs = 120)
+mypairs <- eye.pairs.from.fourcorner(450)
+mypairs <- eye.pairs.from.fourcorner(num_pairs = 450)
+mypairs <- eye.pairs.from.fourcorner(num_pairs = 450, mux = 0, muy = 0,
+           sigmax = 1, sigmay = 1, rhoxy = 0) # as for May's matrices
 }
 
diff --git a/man/eye.pairs.from.normal.Rd b/man/eye.pairs.from.normal.Rd
index 7e317ce..f2b1ea2 100644
--- a/man/eye.pairs.from.normal.Rd
+++ b/man/eye.pairs.from.normal.Rd
@@ -3,11 +3,11 @@
 \alias{eye.pairs.from.normal}
 \title{Sample pairs of interactions from a bivariate normal distribution}
 \usage{
-eye.pairs.from.normal(NumPairs = 10, mux = -1, muy = 0.5, sigmax = 1/4,
+eye.pairs.from.normal(num_pairs = 10, mux = -1, muy = 0.5, sigmax = 1/4,
   sigmay = 1/4, rhoxy = -2/3)
 }
 \arguments{
-\item{NumPairs}{The number of pairs to sample}
+\item{num_pairs}{The number of pairs to sample}
 
 \item{mux}{The desired mean for the first marginal distribution}
 
@@ -17,16 +17,19 @@ eye.pairs.from.normal(NumPairs = 10, mux = -1, muy = 0.5, sigmax = 1/4,
 
 \item{sigmay}{The desired standard deviation for the second marginal distribution}
 
-\item{rhoxy}{The desired correlation}
+\item{rhoxy}{The desired correlation between the pairs of interaction strengths}
 }
 \value{
-A NumPairs x 2 matrix of interaction strengths sampled from the distribution
+A num_pairs x 2 matrix of interaction strengths sampled from the distribution
 }
 \description{
 Sample pairs of interactions from a bivariate normal distribution
 }
 \examples{
 mypairs <- eye.pairs.from.normal() # default values
-mypairs <- eye.pairs.from.normal(NumPairs = 120)
+mypairs <- eye.pairs.from.normal(450)
+mypairs <- eye.pairs.from.normal(num_pairs = 450)
+mypairs <- eye.pairs.from.normal(num_pairs = 450, mux = 0, muy = 0,
+           sigmax = 1, sigmay = 1, rhoxy = 0) # as for May's matrices
 }
 
diff --git a/man/eye.parameterize.M.Rd b/man/eye.parameterize.M.Rd
new file mode 100644
index 0000000..ce6e1a4
--- /dev/null
+++ b/man/eye.parameterize.M.Rd
@@ -0,0 +1,40 @@
+% Generated by roxygen2 (4.0.1): do not edit by hand
+\name{eye.parameterize.M}
+\alias{eye.parameterize.M}
+\title{Builds a community matrix M based on a food web, and a parameterization for the bivariate distribution}
+\usage{
+eye.parameterize.M(FW, distribution_pairs = "Normal", mux = -1, muy = 0.5,
+  sigmax = 1/4, sigmay = 1/4, rhoxy = -2/3, mu_diagonal = 0,
+  sigma_diagonal = 0)
+}
+\arguments{
+\item{FW}{A food web object, produced by the method eye.foodweb.cascade or similar}
+
+\item{distribution_pairs}{A bivariate distribution for the pairs (possibilities are "Normal", "FourCorner", or an empirical distribution, represented by a Nx2 matrix of pairs)}
+
+\item{mux}{Mean of the negative interaction strengths}
+
+\item{muy}{Mean of the positive interaction strengths}
+
+\item{sigmax}{Standard deviation of the negative interaction strengths}
+
+\item{sigmay}{Standard deviation of the positive interaction strengths}
+
+\item{rhoxy}{Correlation between the interaction strengths}
+
+\item{mu_diagonal}{(0 by default) mean of the diagonal elements}
+
+\item{sigma_diagonal}{(0 by default) standard deviation for the diagonal elements}
+}
+\value{
+The community matrix M
+}
+\description{
+Builds a community matrix M based on a food web, and a parameterization for the bivariate distribution
+}
+\examples{
+M <- eye.parameterize.M(eye.foodweb.cascade())
+M <- eye.parameterize.M(eye.foodweb.niche(size = 30, connectance = 0.2),
+     distribution_pairs = "FourCorner", mux = -1, muy = 1, sigmax = 2, sigmay = 1, rhoxy = -0.5)
+}
+
diff --git a/man/eye.plot.approximations.Rd b/man/eye.plot.approximations.Rd
new file mode 100644
index 0000000..6e46545
--- /dev/null
+++ b/man/eye.plot.approximations.Rd
@@ -0,0 +1,24 @@
+% Generated by roxygen2 (4.0.1): do not edit by hand
+\name{eye.plot.approximations}
+\alias{eye.plot.approximations}
+\title{Plot a graph showing the quality of the approximations}
+\usage{
+eye.plot.approximations(approximation, only_bulk = FALSE)
+}
+\arguments{
+\item{approximation}{An object returned by eye.approximate.ReL1}
+
+\item{only_bulk}{If TRUE, plot only the eigenvalues close to the center of the distribution}
+}
+\value{
+A ggplot plot
+}
+\description{
+Plot a graph showing the quality of the approximations
+}
+\examples{
+M <- eye.parameterize.M(eye.foodweb.cascade())
+Appr <- eye.approximate.ReL1(M)
+eye.plot.approximations(Appr)
+}
+