[z_Deprecated] Evaluation

Comprehensive evaluation is implemented in the script /benchmark/evaluate.py, with partial visualizations for the real-world datasets implemented in /benchmark/evaluate_alt_rw_plots.py.

Similarly to the methods benchmarking, for simulation, it can be wrapped and run in command line using options 2, 3 in /benchmark/evaluate.py with the iterative_methods_running_evaluate function with run_or_evaluate = 'evaluate'. On the other hand, for real-world dataset, it can be wrapped and run in command line using options 4, 5 in /benchmark/evaluate.py with the iterative_methods_running_evaluate function with run_or_evaluate = 'evaluate'. For reproducibility, the line plots for real-world evaluations are implemented in /benchmark/evaluate_alt_rw_plots.py.

The evaluation section is divided into two sections: 3.1 generating summary statistics for each simulated dataset or for each real-world dataset; and 3.2 generating intuitive visualizations based on the summary statistics from 3.1.

3.1 Generating summary statistics for each simulated dataset or for each real-world dataset

3.1.1 Simulated datasets

To aid your intuitive understanding, you can do the following for a minimally viable run in iPython. This will generate a folder of files similar to benchmark/trial/simulation_data_eval_count_norelation/out_1.25_0.5_0.5_iter_1, which has output subfolders consisting of evaluation results from each method, a csv file with global summary stats, and a pdf figure for visualizing cross-study bias/batch effect removal across all methods.

• In ipython:

from evaluate import *
overall_path = './trial' # pls make sure this dir exists
or_l = [1.25]
cond_effect_val_l = [0.5]
batch_effect_val_l = [0.5]
GLOBAL_DATATYPE = 'count'
related = 'no'
binarizing_agent_biovar = 'cond_1'
iterative_methods_running_evaluate(run_or_evaluate = 'evaluate', datatype = GLOBAL_DATATYPE, iter = 1, or_l = or_l, cond_effect_val_l = cond_effect_val_l, batch_effect_val_l = batch_effect_val_l, 
                    address_XY_dir_path = './trial/simulate', 
                    output_dir_path = f'./trial/simulation_data_output_{GLOBAL_DATATYPE}_{related}relation', 
                    eval_dir_path = overall_path+f"/simulation_data_eval_{GLOBAL_DATATYPE}_{related}relation",
                    methods_list = methods_list_dict[GLOBAL_DATATYPE], 
                    binarizing_agent_biovar = binarizing_agent_biovar)

This will generate the evaluation files, as exemplified in the folder here. It should contain one folder per benchmarked method, a summary stats csv file named 'global_benchmarking_stats....csv' and a PCoA visualization in pdf format.

Alternatively, we provide a command line option, which might require you to go in and change certain codes in evaluate.py, but allows you a more streamlined running experience:

1. Replace the parameter combination set up following the line ### parameter combination set-up in evaluate.py.
2. Ensure that the methods run and the paths where the input files are taken and output files generated are correct following the line elif args['option'] == 2: in evaluate.py.
3. Run python3 /benchmark/evaluate.py -o 2 -d relab -r yes -i 1 -a cond_1 -p /athena/linglab/scratch/chf4012. You can easily replace relab with count (the default) or yes with no (default), 1 with other iteration numbers, and the path with the overall path where all your files are generated. For details, please check argparse of the file evaluate.py with python3 evaluate.py -h.

For this alternative option, you can scale run with slurm by referencing the file /benchmark/slurm_bash_scripts/step2_evaluate_run_sim.sh, revising it, and moving it to benchmark folder to slurm run.

3.1.2 Real-world dataset

For the iPython option, you can do the following for a minimally viable run in iPython. This will generate a folder of files similar to benchmark/trial/autism_2_microbiomeHD_eval, which has output subfolders consisting of evaluation results from each method, a csv file with global summary stats, and a pdf figure for visualizing batch effect removal across all methods.

• In ipython

# autism 2 microbiomeHD
output_dir_path = './trial/autism_2_microbiomeHD_eval'
address_directory = '../data/cleaned_data/autism_2_microbiomeHD'
vars_use = ["Dataset"]
IDCol = 'Sam_id'

address_Y = "../data/cleaned_data/autism_2_microbiomeHD/autism_2_microbiomeHD_meta_data.csv"

# nobc
address_X = "../data/cleaned_data/autism_2_microbiomeHD/autism_2_microbiomeHD_count_data.csv"
data_mat, meta_data = load_results_from_benchmarked_methods(address_X, address_Y)
### nobc
Evaluate(data_mat, meta_data, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_nobc/autism_2_microbiomeHD_nobc', "DiseaseState", 30, [], 'Sam_id', method = 'nobc')

### harmony
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_harmony_adjusted_count.csv"
res_h, meta_data_h = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(res_h, meta_data_h, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_harmony/autism_2_microbiomeHD_harmony', "DiseaseState", 30, [], 'Sam_id', method = 'harmony')

# benchmarking other methods: 
### combat (combat_seq)
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_combat_seq.csv"
data_mat_combat, meta_data_combat = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(data_mat_combat, meta_data_combat, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_combat_seq/autism_2_microbiomeHD_combat_seq', "DiseaseState", 30, [], 'Sam_id', method = 'combat_seq')

### limma
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_limma.csv"
data_mat_limma, meta_data_limma = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(data_mat_limma, meta_data_limma, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_limma/autism_2_microbiomeHD_limma', "DiseaseState", 30, [], 'Sam_id', method = 'limma')

### MMUPHin
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_MMUPHin.csv"
data_mat_mmuphin, meta_data_mmuphin = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(data_mat_mmuphin, meta_data_mmuphin, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_MMUPHin/autism_2_microbiomeHD_MMUPHin', "DiseaseState", 30, [], 'Sam_id', method = 'MMUPHin')

### ConQuR
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_ConQuR.csv"
data_mat_conqur, meta_data_conqur = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(data_mat_conqur, meta_data_conqur, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_ConQuR/autism_2_microbiomeHD_ConQuR', "DiseaseState", 30, [], 'Sam_id', method = 'ConQuR')

### ConQuR_libsize
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_ConQuR_libsize.csv"
data_mat_conqur_libsize, meta_data_conqur_libsize = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(data_mat_conqur_libsize, meta_data_conqur_libsize, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_ConQuR_libsize/autism_2_microbiomeHD_ConQuR_libsize', "DiseaseState", 30, [], 'Sam_id', method = 'ConQuR_libsize')

### percentile_norm
address_X = "./trial/rw_methods_out/autism_2_microbiomeHD/autism_2_microbiomeHD_percentile_norm.csv"
data_mat_percentile_norm, meta_data_percentile_norm = load_results_from_benchmarked_methods(address_X, address_Y)
Evaluate(data_mat_percentile_norm, meta_data_percentile_norm, 'Dataset', output_dir_path + '/output_autism_2_microbiomeHD_percentile_norm/autism_2_microbiomeHD_percentile_norm', "DiseaseState", 30, [], 'Sam_id', method = 'percentile_norm')

### global evaluation
input_frame_path = "../data/cleaned_data/autism_2_microbiomeHD/autism_2_microbiomeHD_count_data.csv"
bio_var = "DiseaseState"
dataset_name = "autism_2_microbiomeHD"
methods_list = ["nobc", "harmony", "combat_seq", "limma", "MMUPHin", "ConQuR", "ConQuR_libsize", "percentile_norm"]
global_eval_dataframe(input_frame_path, bio_var, dataset_name, methods_list, "./trial/rw_methods_out", './trial/autism_2_microbiomeHD_eval', simulate = False)

### multi-method plot
df_l = [data_mat, res_h, data_mat_combat, data_mat_limma, data_mat_mmuphin, data_mat_conqur, data_mat_conqur_libsize, data_mat_percentile_norm]
methods = ["nobc", "harmony", "combat_seq", "limma", "MMUPhin", "ConQuR", "ConQuR_libsize", "percentile_norm"]
meta_data_l = [meta_data, meta_data_h, meta_data_combat, meta_data_limma, meta_data_mmuphin, meta_data_conqur, meta_data_conqur_libsize, meta_data_percentile_norm]
plot_PCOA_multiple('autism_2_microbiomeHD', df_l, methods, meta_data_l, used_var="Dataset", output_root= output_dir_path + '/')

3.2 Generating intuitive visualizations based on the summary statistics

After generating summary statistics over 1000 iterations in simulations, or across all real-world datasets, it is important to visualize in an intuitive manner to allow the users to understand what methods are better in what scenarios. Different from the other sections, this section requires users to either have run all iterations of simulations, or to have run all the real-world datasets of interest. This is because these summary plots would not make sense without all the iterations/diversity of real-world datasets.

This section is made easy with the option 3 of evaluate.py for simulation dataset and directly running python3 evaluate_alt_rw_plots and for real-world datasets respectively. These generate line plots to visualize the statistics of different batch corrected methods in various simulation settings and in various real-world datasets. If you seek to run these codes, please double check the part of evaluate.py following line elif args['option'] == 3: and lines 299-311 of evaluate_alt_rw_plots carefully to ensure the path of the files are correct. Note that there is a certain level of customness to the codes to ensure the generated figures are good-looking.

For example, one can run the following code for simulation:

python3 evaluate.py -o 3 -d relab -r yes

This should reproduce our plots as exemplified in the folder here, one set of figures for each parameter setting.

And the following code for the real-world data:

python3 evaluate_alt_rw_plots.py

This should generate plots as exemplified ini this folder.