Mean RGB image out of 4D np.array
def average_image(a_lot_of_images):
    avg = np.zeros((a_lot_of_images.shape[1], a_lot_of_images.shape[2], a_lot_of_images.shape[3]), dtype=float)
    for i in range(a_lot_of_images.shape[0]):
        avg[:, :, 0] += a_lot_of_images[i, :, :, 0]
        avg[:, :, 1] += a_lot_of_images[i, :, :, 1]
        avg[:, :, 2] += a_lot_of_images[i, :, :, 2]
    for i in [0,1,2]:
        avg[:, :, i] = avg[:, :, i]/a_lot_of_images.shape[0]
    avg = np.array(np.round(avg), dtype=np.uint8)
    return avg
works! Is there a way to do this in a more elegant way?
I have a 4D data np.array consisting of 150 images 50X70 decomposed to 3 channels.
The shape of the data is (150, 50, 70, 3). 
I need the average of those 150 images by channel, with the shape (50, 70, 3) 
(where R channel of average would be the average of 150 R channels, etc.)
I've tried:
average = data.mean(axis=0)
averageimage = Image.fromarray(average, 'RGB')
averageimage.show()
but even tho it gives correct shape, the image looks like random color noise. 
EDIT: I've tried 
def average_image(a_lot_of_images):
    avg = np.zeros((a_lot_of_images.shape[1], a_lot_of_images.shape[2], a_lot_of_images.shape[3]))
    for i in range(a_lot_of_images.shape[0]):
        avg[:,:,0] += a_lot_of_images[i,:,:,0]
        avg[:,:,1] += a_lot_of_images[i,:,:,1]
        avg[:,:,2] += a_lot_of_images[i,:,:,2]
    for i in [0,1,2]:
        avg[:,:,i] = avg[:,:,i]/a_lot_of_images.shape[0]
    return avg
and the output still looks like color noise.