保存したVAEの学習済みモデルを新しいデータ（normal画像、anomal画像）に使いたいと考えています。
以下の通りload_modelで読み込んだのですが、「Cannot create group in read only mode.」というエラーが出てうまく読み込む事が出来ませんでした。

色々ネットで調べたところ、モデルのcompile前に読み込みを行うとこのようなエラーが出てしまうようでした。
ですが、以下のコードでは読み込み前にcompileを行っており、どうすればエラーが解消されるのかが分かりません。
どなたかアドバイスして頂けないでしょうか。
お手数をおかけしますがどうぞよろしくお願い致します。

from __future__ import absolute_import
from __future__ import division
from __future__ import print_function

from keras.layers import Lambda, Input, Dense
from keras.models import Model
from keras.models import Sequential, model_from_json
from keras.losses import mse, binary_crossentropy
from keras.layers import Conv2D, Flatten, Lambda
from keras.layers import Reshape, Conv2DTranspose
from keras.utils import plot_model, np_utils 
from keras.utils import plot_model
from keras import optimizers
from keras import backend as K
from keras.preprocessing.image import array_to_img, img_to_array,load_img
from keras.preprocessing.image import ImageDataGenerator
import numpy as np
import matplotlib.pyplot as plt
import argparse
import os
import re
import glob
import random as rn
import tensorflow as tf
import cv2
from PIL import Image
from google.colab.patches import cv2_imshow
from keras.models import load_model
import numpy as np
from sklearn import metrics
from sklearn.metrics import mean_squared_error
from sklearn.metrics import roc_curve
from sklearn.metrics import roc_auc_score
import h5py

import warnings
warnings.filterwarnings('ignore')

%matplotlib inline
get_ipython().run_line_magic('matplotlib', 'inline')

def sampling(args):
    """Reparameterization trick by sampling fr an isotropic unit Gaussian.

    # Arguments
        args (tensor): mean and log of variance of Q(z|X)

    # Returns
        z (tensor): sampled latent vector
    """

    z_mean, z_log_var = args
    batch = K.shape(z_mean)[0]
    dim = K.int_shape(z_mean)[1]
    # by default, random_normal has mean=0 and std=1.0
    epsilon = K.random_normal(shape=(batch, dim))
    return z_mean + K.exp(0.5 * z_log_var) * epsilon

def plot_results(models,
                 data,
                 batch_size=32,
                 model_name="vae_OCT"):
    """Plots labels and MNIST digits as function of 2-dim latent vector

    # Arguments
        models (tuple): encoder and decoder models
        data (tuple): test data and label
        batch_size (int): prediction batch size
        model_name (string): which model is using this function
    """

    encoder, decoder = models
    x_test = data
    os.makedirs(model_name, exist_ok=True)

    filename = os.path.join(model_name, "vae_mean.png")

#normal画像
filenames = glob.glob("./right_ver_normal_resize/*.tif")
test_normal= []

for filename in filenames:
    img = img_to_array(load_img(
    filename, color_mode = "grayscale"
    , target_size=(512,496)))
    test_normal.append(img)

test_normal = np.asarray(test_normal )

#anomaly画像
filenames = glob.glob("./red_by_ueno_resize/*.tif")
test_anomaly= []

for filename in filenames:
    img = img_to_array(load_img(
    filename, color_mode = "grayscale"
    , target_size=(512,496)))
    test_anomaly.append(img)

test_anomaly = np.asarray(test_anomaly)

image_size = (512, 496)# X.shape[1]
original_dim = 512 * 496
test_normal = np.reshape(test_normal, [-1, original_dim])# test_normal = np.reshape(test_normal, [-1, image_size_width,image_size_height,1])
test_anomaly = np.reshape(test_anomaly, [-1, original_dim])# test_anomaly = np.reshape(test_anomaly, [-1, image_size_width,image_size_height,1])
test_normal = test_normal.astype('float32') / 255
test_anomaly= test_anomaly.astype('float32') / 255

print(test_normal.shape)
print(test_anomaly.shape)

# network parameters
input_shape = (original_dim,)
batch_size = 50
intermediate_dim = 128
latent_dim = 2
epochs = 5

# VAE model = encoder + decoder
# build encoder model
inputs = Input(shape=input_shape, name='encoder_input')
x = Dense(intermediate_dim, activation='relu')(inputs)
z_mean = Dense(latent_dim, name='z_mean')(x)
z_log_var = Dense(latent_dim, name='z_log_var')(x)

# use reparameterization trick to push the sampling out as input
# note that "output_shape" isn't necessary with the TensorFlow backend
z = Lambda(sampling, output_shape=(latent_dim,), name='z')([z_mean, z_log_var])

# instantiate encoder model
encoder = Model(inputs, [z_mean, z_log_var, z], name='encoder')
encoder.summary()
plot_model(encoder, to_file='vae_mlp_encoder.png', show_shapes=True)

# build decoder model
latent_inputs = Input(shape=(latent_dim,), name='z_sampling')
x = Dense(intermediate_dim, activation='relu')(latent_inputs)
outputs = Dense(original_dim, activation='sigmoid')(x)

# instantiate decoder model
decoder = Model(latent_inputs, outputs, name='decoder')
decoder.summary()
plot_model(decoder, to_file='vae_mlp_decoder.png', show_shapes=True)

# instantiate VAE model
outputs = decoder(encoder(inputs)[2])
vae = Model(inputs, outputs, name='vae_mlp')


if __name__ == '__main__':
    parser = argparse.ArgumentParser()
    help_ = "Load h5 model trained weights"
    parser.add_argument("-w", "--weights", help=help_)
    help_ = "Use mse loss instead of binary cross entropy (default)"
    parser.add_argument("-m", "--mse", help=help_, action='store_true')
    args = parser.parse_args([])# args =  ["--weights", "--mse"]
    models = (encoder, decoder)
    #data = (x_test)



    # VAE loss = mse_loss or xent_loss + kl_loss
    if args.mse:
        reconstruction_loss = mse(K.flatten(inputs), K.flatten(outputs))
    else:
        reconstruction_loss = binary_crossentropy(K.flatten(inputs),
                                                  K.flatten(outputs))

    reconstruction_loss *= original_dim
    kl_loss = 1 + z_log_var - K.square(z_mean) - K.exp(z_log_var)
    kl_loss = K.sum(kl_loss, axis=-1)
    kl_loss *= -0.5
    vae_loss = K.mean(reconstruction_loss + kl_loss)
    vae.add_loss(vae_loss)
    Adam = optimizers.Adam(lr=0.0005)
    vae.compile(optimizer=Adam)
    vae.summary()
    plot_model(vae, to_file='vae_cnn.png', show_shapes=True)

# load model
vae = load_model("model.ep03.h5")

test_normal_pred = vae.predict(test_normal, verbose=1)
test_anomaly_pred = vae.predict(test_anomaly, verbose=1)

score = []
for i in range(len(test_normal)):
  squared_error_normal = mean_squared_error(test_normal[i,:], test_normal_pred[i,:])

for i in range(len(test_anomaly)):
  squared_error_anomaly = mean_squared_error(test_anomaly[i,:], test_anomaly_pred[i,:])
  score.append(squared_error_anomaly)

score = np.array(score)
#print(score.shape)
#print(score)

#ROC曲線の描画
y_true = np.zeros(len(test_normal)+len(test_anomaly))
y_true[len(test_normal):] = 1#0:正常、1：異常
# FPR, TPR(, しきい値) を算出
fpr, tpr, thresholds = roc_curve(y_true, score)

# AUC
auc = roc_auc_score(y_true, score)

# ROC曲線をプロット
plt.figure
plt.plot(fpr, tpr, label='VAE method (area = %.2f)'%auc, c="r")
plt.legend()
plt.title('ROC curve')
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.grid(True)
plt.savefig("./ROC curve.png")
plt.show()