Is there a way to automatically make a "collage" of plots with matplotlib?

ax.set_adjustable('box-forced')
ax2.set_adjustable('box-forced')

import numpy as np

def unlink_wrap(dat, lims=[-np.pi, np.pi], thresh = 0.95):
    """
    Iterate over contiguous regions of `dat` (i.e. where it does not
    jump from near one limit to the other).

    This function returns an iterator object that yields slice
    objects, which index the contiguous portions of `dat`.

    This function implicitly assumes that all points in `dat` fall
    within `lims`.

    """
    jump = np.nonzero(np.abs(np.diff(dat)) > ((lims[1] - lims[0]) * thresh))[0]
    lasti = 0
    for ind in jump:
        yield slice(lasti, ind + 1)
        lasti = ind + 1
    yield slice(lasti, len(dat))

x = np.arange(0, 100, .1)
y = x.copy()

lims = [0, 24]

x = (x % lims[1])

fig, ax = matplotlib.pyplot.subplots()

for slc in unlink_wrap(x, lims):
    ax.plot(x[slc], y[slc], 'b-', linewidth=2)

ax.plot(x, y, 'r-', zorder=-10)
ax.set_xlim(lims)

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.cm as mplcm
import matplotlib.colors as colors
import numpy as np
from random import randint

fig = plt.figure()
ax = Axes3D(fig)

num_lines = 50
lines, distances = [], []
for x in range(0, num_lines):
    xs = [0, randint(10, 150)]
    ys = [0, 5]
    zs = [0, randint(10,60)]
    lines.append([xs, ys, zs])

    distance = np.sqrt((xs[0]-xs[1])**2+(ys[0]-ys[1])**2+(zs[0]-zs[1])**2)  # calculate length of a line
    distances.append(distance)

distances_normalised = distances / (np.max(distances)/num_lines)  # normalise lengths of lines
cNorm = colors.Normalize(vmin=0, vmax=num_lines - 1)

cm = plt.get_cmap('hot_r')
scalarMap = mplcm.ScalarMappable(norm=cNorm, cmap=cm)

for dist, line in enumerate(lines):
    ax.plot(line[0], line[1], line[2], color=scalarMap.to_rgba(distances_normalised[dist]))

ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')

plt.show()

f, axarr = plt.subplots(nrows = 1, ncols=3)

In [7]: axarr.shape
Out[7]: (3,)

In [4]: axarr
Out[4]: 
array([<matplotlib.axes._subplots.AxesSubplot object at 0x7f0e299dead0>,
       <matplotlib.axes._subplots.AxesSubplot object at 0x7f0e107d3110>,
       <matplotlib.axes._subplots.AxesSubplot object at 0x7f0e1804d5d0>], dtype=object)

In [5]: axarr[0,0]
---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-5-3e75475f2c86> in <module>()
----> 1 axarr[0,0]

IndexError: too many indices

proc_plot <- function(sub) {

  data_Group1 <- sub[sub$Groups == "Group1", ]
  data_Group2 <- sub[sub$Groups == "Group2", ]

  min_rdg <- min(data_Group1$Readings, data_Group2$Readings)
  max_rdg <- max(data_Group1$Readings, data_Group2$Readings)

  # Group1
  graph_Group1 <- ggplot(data_Group1, aes(x = TimePoint, y = Readings, group = Scale)) +
    labs(title="Group1", x="Time point", y="Readings") +
    scale_x_continuous(breaks = c(42.5, 47.5, 52.5, 57.5, 62.5, 67.5, 72.5),
                       labels = c("1", "2", "3", "4", "5", "6", "7")) +
    geom_line(aes(color = Scale), na.rm = TRUE) +
    geom_point(aes(color = Scale),size = 1.5, na.rm = TRUE) +
    scale_color_continuous(name = "Scale", breaks = c(5, 10, 25, 50, 75, 90)) +
    scale_y_continuous(limits = c(min_rdg, max_rdg)) +
    theme(legend.key.height = unit(3.5, "cm"))

  # Group2
  graph_Group2 <- ggplot(data_Group2, aes(x = TimePoint, y = Readings, group = Scale)) +
    labs(title="Group2", x="Time point", y="Readings") +
    scale_x_continuous(breaks = c(42.5, 47.5, 52.5, 57.5, 62.5, 67.5, 72.5),
                       labels = c("1", "2", "3", "4", "5", "6", "7")) +
    geom_line(aes(color = Scale), na.rm = TRUE) +
    geom_point(aes(color = Scale), size = 1.5, na.rm = TRUE) +
    scale_color_continuous(name = "Scale", breaks = c(5, 10, 25, 50, 75, 90)) +
    scale_y_continuous(limits = c(min_rdg, max_rdg)) +
    theme(legend.key.height = unit(3.5, "cm"))

  png (paste0("ScatterPlot_", sub$Biomarkers[[1]], ".png"), height=600, width=1111)
    output <- grid.arrange(graph_Group1, graph_Group2, nrow = 1, 
                           top=textGrob(sub$Biomarkers[[1]], gp=gpar(fontsize=20)))
  dev.off()

  return(output)
}

# BUILD PLOT LIST AND PNG FILES
plot_list <- by(data, data$Biomarkers, proc_plot)

dev.off()
grid.draw(plot_list$Gene1)

dev.off()
grid.draw(plot_list$Gene2)