is this code crrectly disply spearman matrix for sensitivity analysis?

import os
import shutil
import pickle
import time
import numpy as np
import scipy.stats as stats
import matplotlib.pyplot as plt
from smt.sampling_methods import LHS

from scipy.optimize import lsq_linear
from collections import OrderedDict
from functions_py3 import *

# Total time - start
total_time_start = time.time()


# Measured frequencies to perform FE model updating on
mode_list = ['Mode 01', 'Mode 02', 'Mode 03', 'Mode 04', 'Mode 05', 'Mode 06','Mode 07']


# -------------------------------------------------------------------------- #
# INITIAL
# -------------------------------------------------------------------------- #

# Folder structure
folder_1 = '02_Analysis files'
folder_2 = '03_Results/01_Perturbed'
folder_3 = '04_Figures'

folders = [folder_1, folder_2, folder_3]

# Establish folder structure
for folder in folders:

    # Remove existing folders with content
    if os.path.isdir(folder):
        shutil.rmtree(folder)
        time.sleep(1)

    # Create folder
    os.makedirs(folder)

# Import measured frequencies and mode shapes
freqs_est = np.loadtxt('01_Initial/02_Measured/NB_0_frequencies_all.txt')
modes_est = np.loadtxt('01_Initial/02_Measured/NB_0_modes_all.txt')

# Copy relevant model files to directory
src = ('01_Initial/01_FE model/')

files_to_copy = ['NB.cae', 'NB.jnl']

# Copy
for file in files_to_copy:

    shutil.copy(src + file, file)

# Set work directory
cwd = os.getcwd()
os.chdir(cwd)

# ---------------------------------------------------------------------------#
# ABAQUS MODEL
# ---------------------------------------------------------------------------#

n_samples = 800  # Number of samples


LB = np.array([0.9*7850, 0.8*2.1E11, 0.8*2.1E11, 0.8*2.1E11, 0.8*2.1E11, 0.8*2.1E11, 
               0.8*2.1E11, 0.9*7850,  0.8*2.1E11, 0.1*2E7, 0.1*1E9, 0.1*2E8, 0.1*5E6, 0.5*40])




UB = np.array([1.1*7850, 1.2*2.1E11, 1.2*2.1E11, 1.2*2.1E11, 1.2*2.1E11, 1.2*2.1E11, 
               1.2*2.1E11, 1.1*7850,  1.2*2.1E11, 100*2E7, 100*1E9, 100*2E8, 100*5E6, 2*40])




xlimits = np.vstack((LB, UB)).T

sampling = LHS(xlimits=xlimits, criterion='ese')
LHS_matrix = sampling(n_samples)

Y = np.zeros((LHS_matrix.shape[0], 2*len(mode_list)))

for j in range(n_samples):
    theta_0 = LHS_matrix[j,:]    
    
    # Job name
    job_name = 'Analysis1'
    with open('job_name_list' + '.txt', 'w') as fout:
        fout.write(job_name)
    
    # Parameter list
    np.save('parameter_list' + '.npy', theta_0, allow_pickle=True,
            fix_imports=True)
    
    # ---------------- #
    # RUN ANALYSIS
    # ---------------- #
    
    # Variables
    script_name = 'abaqus_NB_upd'
    
    # Run script
    print('--------------------------------------------')
    print('Initial analysis ({})'.format(j+1))
    t0 = time.time()
    
    os.system('abaqus cae noGUI=' + script_name)
    
    t1 = time.time()
    print('Initial analysis ({}) done - '.format(j+1) + str(round(t1 - t0, 3)) + ' sec.')
    print('--------------------------------------------')
    print()
    
    # Wait 3 sec for the Abaqus analysis to properly close
    time.sleep(3)
    
    # ------------------------------------ #
    # POST-PROCESS
    # ------------------------------------ #
    
    # Import
    import_folder = '03_Results/'
    
    freqs_num = np.load(import_folder + job_name + '_frequencies_all.npy')
    modes_num = np.load(import_folder + job_name + '_modes_all.npy')
    
    # Sorted MAC and MMI
    (MAC_initial_model, MMI_initial_model, results_MAC, results_MMI) = get_MAC_MMI(freqs_est, modes_est,
                                                         freqs_num, modes_num)
    
    # General parameters
    rows_est, cols_est = np.shape(modes_est)
    rows_num, cols_num = np.shape(modes_num)
    
    
    # ---------------------------------------------------------------------------#
    # FEM UPDATING
    # ---------------------------------------------------------------------------#
    
    # ------------------------------------ #
    # INITIALIZATION
    # ------------------------------------ #
    
    # Variables
    q = len(mode_list)
    
    # Declaration of variables
    # Variables saved for each iteration
    lamda = np.zeros((q), dtype=float)
    mac = np.zeros((q), dtype=float)
    lamda_all = np.zeros((len(freqs_est)), dtype=float)
    mac_all = np.zeros((len(freqs_est)), dtype=float)
    lamda_m = np.zeros(len(mode_list), dtype=float)
    
    # Start values for the iterations
    for idx, mode in enumerate(mode_list):
    
        # Establish frequencies based on the mode_list string
        string_idx = mode[-2:]
    
        # Correct for 0 in the string
        if string_idx[0] == 0:
    
            f_idx = int(mode[-1]) - 1
    
        else:
            f_idx = int(string_idx) - 1
    
        # Measured frequencies (to perform FEM updating on)
        lamda_m[idx] = freqs_est[:, 1][f_idx]
    
        # Numerical frequencies (to perform FEM updating on)
        # Top match based on MMI
        # lamda[0, idx] = results_MAC[mode][0, 2]
        lamda[idx] = results_MMI[mode][0, 2]
        mac[idx] = results_MMI[mode][0, 1]
    
    # All numerical frequencies
    for idx, mode in enumerate(results_MMI.keys()):
        
        # lamda_all[0, idx] = results_MAC[mode][0, 2]
        lamda_all[idx] = results_MMI[mode][0, 2]
        mac_all[idx] = results_MMI[mode][0, 1] 
    
    Y[j,:len(lamda)] = lamda
    Y[j, len(lamda):len(lamda)+len(mac)] = mac
    
    # Wait 5 sec for the Abaqus analysis to properly close
    # print('Sleeping 5 sec for Abaqus to properly close.')
    # time.sleep(5)


# The samples with a MAC value lower than 0.50 have been excluded from the total amount of samples.
# Additionally, the samples with a correlation coefficient within the range −0.30 to +0.30 were removed.
# middle_index = Y.shape[1] // 2
# mask = (Y[:, middle_index:] >= 0.5).all(axis=1)
# Y_mod = Y[mask]
LHS_mod = LHS_matrix
Rho = np.zeros((len(LB), 2*len(mode_list)))
for i in range(len(LB)):
    for j in range(2*len(mode_list)):
        Rho[i, j], _ = stats.spearmanr(LHS_mod[:, i], Y_mod[:, j])



rho= Rho.T

# columns = np.any(np.abs(rho) > 0.3, axis=0)
# print(f'Columns with any absolute value larger than 0.3: {np.where(columns)[0]}')

# Plot the heat figure
fig, ax = plt.subplots()  # You can modify the values as per your requirement
# idx = (rho < -0.30) | (rho > 0.30)  # Find the indices of Rho outside the range
cax = ax.imshow(rho,cmap='coolwarm', vmin=-1, vmax=1)  # Only show the values outside the range
ax.set_xticks(np.arange(len(LB)))
ax.set_yticks(np.arange(2*len(mode_list)))
ax.set_xticklabels([r'$\rho_{st}$', r'$E_b$', r'$E_c$', r'$E_s$', r'$E_t$', r'$E_f$', r'$E_{lb}$',
                    r'$\rho_r$', r'$E_r$', r'$K_{x}$', r'$K_{y}$', r'$K_{z}$', r'$K_{\phi}$', r'$m_{bd}$'], rotation=45, fontsize=10)

ax.set_yticklabels([r'$f_1$', r'$f_2$', r'$f_3$', r'$f_4$', r'$f_5$', r'$f_6$', r'$f_7$',
                    r'$MAC_1$', r'$MAC_2$', r'$MAC_3$', r'$MAC_4$', r'$MAC_5$', r'$MAC_6$', r'$MAC_7$'])
ax.set_xlabel('Inputs')
ax.set_ylabel('Outputs')
ax.set_title('Spearman Correlation Coefficient matrix')
cbar = fig.colorbar(cax)
plt.tight_layout()
plt.show()