KerasのCNNにおける複数次元の入力方法について

##前提
KerasのCNNで画像を学習し、学習した重みを保存するプログラムを以下のプログラムで作成済みです。プログラムでは"a0"フォルダの画像群と"a1"フォルダの画像群をCNNで学習し、2クラス分類するための重みを保存しています。
a0フォルダの画像ファイル名は[1-a.jpg, 2-a.jpg, ... ,1000-a.jpg]となっており、
a1フォルダの画像ファイル名は[1001-a.jpg, 1002-a.jpg, ... ,2000-a.jpg]となっています。
全てRGB画像ですのでCNNへの入力は3次元になっています。

##実現したいこと
別の画像フォルダ"b0"と"b1"内の画像を"a0"と"a1"の画像とセットにして6次元の入力で学習させたいです。
b0フォルダの画像ファイル名は[1-b.jpg, 2-b.jpg, ... ,1000-b.jpg]となっており、
b1フォルダの画像ファイル名は[1001-b.jpg, 1002-b.jpg, ... ,2000-b.jpg]となっています。
1-a.jpg(3チャンネル)と1-b.jpg(3チャンネル)を1セットとして計6チャンネルで入力→
同様に2-a.jpgと2-b.jpgの6チャンネルを入力→これをフォルダ内全ての画像(1000ペア)で行い学習させるイメージです。

コードが見辛かったら申し訳ありません。
このような6次元の入力は可能なのでしょうか...?
可能でしたら以下のプログラムを改良したコードを回答欄に載せていただけると助かります...

よろしくお願いします。

####学習プログラム

# coding:utf-8                                                                                                             

import keras
from keras.layers import Conv2D, MaxPooling2D, Lambda, Input, Dense, Flatten, BatchNormalization
from keras.utils import np_utils
from keras.models import Sequential
from keras.layers.convolutional import Conv2D, MaxPooling2D
from keras.layers.core import Dense, Dropout, Activation, Flatten
import numpy as np
from sklearn.model_selection import train_test_split
from PIL import Image
import glob
from keras.utils import plot_model
#mpl.use('Agg')                                                                                                            
import matplotlib.pyplot as plt
import tensorflow as tf
from keras import optimizers


folder = ["a0","a1"]
image_size = 224
epoch_size = 100
batch_size = 256

X = []
Y = []
#画像は50×50ピクセルのRGB形式に変換され、それぞれnumpy形式に変換                                                           
for index, name in enumerate(folder):
    dir = "./" + name
    files = glob.glob(dir + "/*.jpg")
    for i, file in enumerate(files):
        image = Image.open(file)
        image = image.convert("RGB")
        image = image.resize((image_size, image_size))
        data = np.asarray(image)
        X.append(data)
        Y.append(index)
#Xは画像データ、Yは正解ラベルのデータ                                                                                      
X = np.array(X)
Y = np.array(Y)
#画像データを0から1の範囲に変換                                                                                            
X = X.astype('float32')
X = X / 255.0

#正解ラベルの形式を変換                                                                                                                                               
Y = np_utils.to_categorical(Y, 2)

# 学習用データとテストデータに分割                                                                                                                                           
#X_train(訓練データ), X_test(テストデータ), y_train(訓練ラベル), y_test(テストラベル)                                                                                                                              
X_train, X_test, y_train, y_test = train_test_split(X, Y, test_size=0.285)

model = Sequential()

model.add(Conv2D(48, 11, strides=(3, 3), activation='relu', padding='same',input_shape=X_train.shape[1:]))
model.add(MaxPooling2D(pool_size=(3, 3), strides=(2, 2)))
model.add(BatchNormalization())
model.add(Conv2D(128, 5, strides=(3, 3), activation='relu', padding='same'))
model.add(MaxPooling2D(pool_size=(3, 3), strides=(2, 2)))
model.add(BatchNormalization())
model.add(Conv2D(192, 3, strides=(1, 1), activation='relu', padding='same'))
model.add(Conv2D(192, 3, strides=(1, 1), activation='relu', padding='same'))
model.add(Conv2D(128, 3, strides=(1, 1), activation='relu', padding='same'))
model.add(MaxPooling2D(pool_size=(3, 3), strides=(2, 2)))
model.add(BatchNormalization())

model.add(Flatten())
model.add(Dense(2048, activation='relu'))
model.add(Dropout(0.5))
model.add(Dense(2048, activation='relu'))
model.add(Dropout(0.5))

model.add(Dense(2, activation='softmax'))

model.summary();
#plot(model, show_shapes=True, show_layer_names=True, to_file='summury.png')                                               

# コンパイル                                                                                                                                                          
model.compile(loss='categorical_crossentropy',optimizer='Adam',metrics=['accuracy'])

#訓練                                                                                                                                                                                   
history = model.fit(X_train, y_train, epochs=epoch_size, batch_size=batch_size, verbose=1, shuffle=True, validation_split=0.25)

#評価 & 評価結果出力                                                                                                                                                                                           
score = model.evaluate(X_test, y_test, verbose=0)
print('Test loss :', score[0])
print('Test accuracy :', score[1])

### save weights                                                                                                           
json_string = model.to_json()
open('alexnet_model.json', 'w').write(json_string)
model.save_weights('alexnet_weights.h5')