咱们是否能够通过气象图来预测降水量呢？明天咱们来应用 CNN 和 LSTM 进行一个乏味的试验。
咱们这里应用荷兰皇家气象研究所 (也称为 KNMI) 提供的凋谢数据集和公共 api，来获取数据集并且构建模型预测当地的降水量。
数据收集
KNMI 提供的数据集，咱们假如气象雷达产生的信号在反射时会被降水（雨、雪、冰雹等）反射。由雷达捕捉的反射信号的强度称为反射率（以 dBZ 计算），咱们能够粗略认为它与该点的降水强度成正比。当通过依据信号强度映射色标将此反射率数据转换为图像时（默认状况下，KNMI 提供的色标为“viridis”，紫色 / 深蓝色示意较低值，黄色示意较高值）产生 图像就像下图显示的那样。咱们每 5 分钟通过 API 以原始格局获取这些数据。然而 API 有一个配额，每小时只能获取 100 张图像。
定义问题
最原始的也是最简略的预测视频中的下一帧的内容的办法是应用 CNN 和 LSTM。咱们是否能够将预测天气雷达的下一个捕捉信号的问题简化为预测视频中的下一帧的问题呢（雷达的讯号也是图像序列）。所以我收集了一些图像序列，并开始试验各种架构的卷积 LSTM 神经网络。每个训练数据点由 36 个间断的雷达原始文件 (对应于距离 5 分钟的 3 小时的测量) 组成。而后将每个数据点分成两局部。前 18 帧用作“特色”(x)，后 18 帧是神经网络在给定前 18 帧的状况下试图预测的内容(y)。在天气预报方面，给定过来 1.5 小时的降水数据，将来 1.5 小时的降水状况(帧距离为 5 分钟，因而 18 帧对应 1.5 小时)。
为什么是卷积 LSTM
如果你对神经网络和深度学习有点相熟，你可能晓得卷积神经网络 (CNN) 在波及剖析或发现图像中的特定特色和形态的工作上体现十分好。而长短期记忆 (LSTM) 神经网络在波及工夫维度 (如工夫序列预测) 和数据序列 (如图像序列、特定工夫范畴内的信号序列等) 的工作上体现十分好。这次要是因为它们有能力学习数据中的长期依赖关系。因而，钻研人员在 2015 年首次提出了一种联合卷积和 LSTM 层的架构，这样能够预测一系列图像中的下一个图像(他们对其进行基准测试的利用之一是降水预测)，所以本文中也是用相似的模型。
数据预处理
咱们应用了近 160 个间断的 36 次雷达扫描序列，咱们应用 h5py 库能够读取并轻松解决原始数据（如从 KNMI 接管的数据是这个格局）并对它们进行预处理。数据点是从 01-01-2019 到当初的随机日期和工夫中筛选的。因为生成的图像的原始尺寸太大，所以将的图像从原始尺寸 (700×765) 放大到 (315×344)。这是模型能够在正当的工夫内训练的最高分辨率，并且在过程中不会有任何的内存溢出问题。而后将每个序列分成两个相等的局部。前 18 帧用作“特色”(x)，后 18 帧是神经网络试图预测的帧(y)(给定前 18 帧)。最初，我将数据集分成两个独自的数据集，别离用于训练(80%) 和验证(20%)。
执行上述所有工作的代码如上面的代码片段所示:
def create_dataset_from_raw(directory_path, resize_to):
resize_width = resize_to[0]
resize_height = resize_to[1]
batch_names = [directory_path + name for name in os.listdir(directory_path) if os.path.isdir(os.path.join(directory_path, name))]
dataset = np.zeros(shape=(len(batch_names),36,resize_height,resize_width)) # (samples, filters, rows = height, cols = width)
for batch_idx,batch in enumerate(batch_names):
files = [x for x in os.listdir(batch) if x != '.DS_Store']
files.sort()
crn_batch = np.zeros(shape=(36, resize_height, resize_width))
for (idx,raster) in enumerate(files):
fn = batch + '/' + raster
img = h5py.File(fn)
original_image = np.array(img["image1"]["image_data"]).astype(float)
img = Image.fromarray(original_image)
# note that here it is (width, heigh) while in the tensor is in (rows = height, cols = width)
img = img.resize(size=(resize_width, resize_height))
original_image = np.array(img)
original_image = original_image / 255.0
crn_batch[idx] = original_image
dataset[batch_idx] = crn_batch
print("Importing batch:" + str(batch_idx+1))
return dataset
def split_data_xy(data):
x = data[:, 0 : 18, :, :]
y = data[:, 18 : 36, :, :]
return x, y
dataset = create_dataset_from_raw('./data/raw/', resize_to=(315,344))
dataset = np.expand_dims(dataset, axis=-1)
dataset_x, dataset_y = split_data_xy(dataset)
X_train, X_val, y_train, y_val = sk.train_test_split(dataset_x,dataset_y,test_size=0.2, random_state = 42)模型
我应用 Tensorflow 和 Keras 框架进行开发。模型基本上是一个自编码器。自编码器是一种神经网络，它试图升高训练数据的维度，对数据进行压缩，而后能够从压缩后潜在空间的散布的近似值中采样，以生成“新”数据。
咱们的模型看起来像这样:
模型共蕴含 9 层 (输出、输入和 7 个暗藏层)。暗藏层在 ConvLSTM2D 层和 BatchNormalization 层之间替换。ConvLSTM2D 层就像简略的 LSTM 层，然而它们的输出和循环转换卷积。ConvLSTM2D 层在保留输出维度的同时，随着工夫的推移执行卷积运算。你能够把它设想成一个简略的卷积层，它的输入被压平，而后作为输出传递到一个简略的 LSTM 层。ConvLSTM2D 层接管模式为(samples, time, channels, rows, cols) 的张量作为输出，输入模式(samples, timesteps, filters, new_rows, new_cols)。所以它们在一段时间内对一系列帧进行运算。
ConvLSTM2D 层之间的 BatchNormalization 层进行归一化操作
对于所有的层(除了输入层)，都应用 LeakyRelu 激活函数，他比 ReLu 好一些，并且和 ReLu 一样快。
该模型采纳二元穿插熵损失函数和 Adadelta 梯度降落优化器进行拟合。因为数据的高维数，Adadelta 会比经典 Adam 优化器有更好的后果。模型训练了 25 个 epoch(之后开始过拟合)。
模型的代码如下所示:
def create_model():
model = Sequential()
model.add(ConvLSTM2D(filters=64, kernel_size=(7, 7),
input_shape=(18,344,315,1),
padding='same',activation=LeakyReLU(alpha=0.01), return_sequences=True))
model.add(BatchNormalization())
model.add(ConvLSTM2D(filters=64, kernel_size=(5, 5),
padding='same',activation=LeakyReLU(alpha=0.01), return_sequences=True))
model.add(BatchNormalization())
model.add(ConvLSTM2D(filters=64, kernel_size=(3, 3),
padding='same',activation=LeakyReLU(alpha=0.01), return_sequences=True))
model.add(BatchNormalization())
model.add(ConvLSTM2D(filters=64, kernel_size=(1, 1),
padding='same',activation=LeakyReLU(alpha=0.01), return_sequences=True))
model.add(Conv3D(filters=1, kernel_size=(3, 3, 3),
activation='sigmoid',
padding='same', data_format='channels_last'))
return model
model = create_model()
model.compile(loss='binary_crossentropy', optimizer='adadelta')
keras.utils.plot_model(model, to_file="model.png", show_dtype=True, show_layer_activations=True, show_shapes=True)
print(model.summary())
epochs = 25
batch_size = 1
#Fit the model
model.fit(
X_train,
y_train,
batch_size=batch_size,
epochs=epochs,
validation_data=(X_val, y_val),
verbose=1,
)后果
在训练模型之后，应用来自验证数据集的示例数据进行测试。模型的输出是 18 个间断的帧(对应于雷达捕捉到的近 1.5 小时的信号)，它返回下一个 18 个预测帧(对应于接下来的 1.5 小时)。
# pick a random index from validation dataset
random_index = np.random.choice(range(len(X_val)), size=1)
test_serie_X = X_val[random_index[0]]
test_serie_Y = y_val[random_index[0]]
first_frames = test_serie_X
original_frames = test_serie_Y
# predict the next 18 fames
new_prediction = model.predict(np.expand_dims(first_frames, axis=0))
new_prediction = np.squeeze(new_prediction, axis=0)
fig, axes = plt.subplots(2, 18, figsize=(20, 4))
# Plot the ground truth frames.
for idx, ax in enumerate(axes[0]):
ax.imshow(np.squeeze(original_frames[idx]), cmap="viridis")
ax.set_title(f"Frame {idx + 18}")
ax.axis("off")
# Plot the predicted frames.
for idx, ax in enumerate(axes[1]):
ax.imshow((new_prediction[idx]).reshape((344,315)), cmap="viridis")
ax.set_title(f"Frame {idx + 18}")
ax.axis("off")
plt.show()将预测的帧与理论状况进行比拟，后果如下图所示。
实在帧与预测帧十分靠近。这种可视化并没有分明地出现降水零碎挪动的工夫维度和方向，因而上面的两个 GIF 动画能够更好地解释模型的输入。
上面是真值的 GIF
上面是预测 GIF
两个动画也十分靠近，模型正确地捕捉到了降水零碎的挪动方向和强度(黄色更强烈，紫色 / 深蓝色强度较低)。
总结
ConvLSTM 将深度学习的两个外围概念联合起来，并取得了很好的成果。本文的残缺我的项目的代码在这里：
https://avoid.overfit.cn/post/da7a8bd23ec14fc6870e1e4b54e85283
作者：Petros Demetrakopoulos