Step 1: Run the REMD simulation Example directory: examples/AMBER_REMD_with_error_analysis/1-simulation-sander Run your REMD simulation. This example is purely for demonstration, for a more complete tutorial on running REMD in AMBER see here. In this example directory I have run a 1 ns REMD simulation on Alanine dipeptide. Step 2a: Extract observables from the trajectory Example directory: examples/AMBER_REMD_with_error_analysis/2-postprocess_ptraj_bootstrap Example script: a-runptraj.sh In AMBER, extracting of observables is typically done by ptraj. The example script extracts dihedral angles from each replica time series using ptraj. Step 2b: Synchronize energy and observable time series Example directory: examples/AMBER_REMD_with_error_analysis/2-postprocess_ptraj_bootstrap Example script: b-synchronize.sh As noted in the manual, Modular reweighting requires that the observable time series and biasing coordinate (in this case, the energy) must be synchronized. In AMBER there are multiple ways to extract out the potential energy. Here we will use the rem.log combined with the ptraj output. Every 10th exchange syncs with a coordinate write. This step is VERY easy to get wrong. Be careful. Step 2c: Bootstrap case resample observable data Example directory: examples/AMBER_REMD_with_error_analysis/2-postprocess_ptraj_bootstrap Example script: c-bootstrap.sh In order to perform error analysis, we will use bootstrap case resampling on the data. Using the bootstrap.py script included in MR, we will create 50 resample sets and (in the following steps) have to analyze each of these resamples. Of course this means that the processing time will be 50 times greater than without resampling. Step 3a: Convert energies to reduced weights Example directory: examples/AMBER_REMD_with_error_analysis/3-MR-conv_reduced_run_WRE Example script: a-convtoreduced.sh Here we will need to convert the potential energies from the rem.log file to a "reduced weight file". This reduced weight file is readable by , "WRE" the WHAM reduced engine. The conversion of time series data to a reduced weight file is what allows Modular Reweighting to be run on ANY equilibrium biased simulation. We will, however, need to run a program to convert the specific biased data (REMD data) to a reduced weight file. This is done here with the convertREMDtoreduced.pl script. Step 3b: Bootstrap resample the reduced weight file Example directory: examples/AMBER_REMD_with_error_analysis/3-MR-conv_reduced_run_WRE Example script: b-bootstrap.sh As with step 2c above, this time series data also needs to be bootstrap resampled. Step 3c: Run WHAM on resamples, convert density of states to Boltzmann distribution Example directory: examples/AMBER_REMD_with_error_analysis/3-MR-conv_reduced_run_WRE Example script: c-runWREon_resamples_convtoBoltzmann.sh Finally we run WHAM (WRE) on the resampled reduced weight files. The output of the WHAM calculation (in this case of REMD) is the energetic density of states. To get a Boltzmann distribution, we simply multiply this by the Boltzmann factor for any temperature of interest. Step 4a: Reweight observables Example directory: examples/AMBER_REMD_with_error_analysis/4-postprocess-MR-probability Example script: a-reweight_observables_using_multihist.sh The density of states and energetic probability distribution (from step 3c) are in-themselves interesting data values, but we often are interested in observables. Here we use the probability data from step 3c to reweight our observable data (resamples) from step 2c using multihist.py. Step 4b: Combine resamples to get distribution with error bars Example directory: examples/AMBER_REMD_with_error_analysis/4-postprocess-MR-probability Example script: b-combine_resamples_with_avg_dev_script.sh Step 4a above gives 50 sets of reweighted distributions, we would like to look at a single distribution with error bars. We can do this using the avg_dev_distribution.py script. ---- Below, optional comparison to raw histograms and plotting --- Step 5a: Extract temperature time series using ptraj Example directory: examples/AMBER_REMD_with_error_analysis/5-compare_to_raw_canonical Example script: a-temperature_time_series.sh It is often worthwhile to check reweighted data against the raw histograms themselves. If want to do this, we start by extracting temperature time series using ptraj. (Since the output trajectories are REPLICA trajectories, not temperature trajectories) Step 5b: Bootstrap resample temperature time series data Example directory: examples/AMBER_REMD_with_error_analysis/5-compare_to_raw_canonical Example script: b-bootstrap.sh To see the error in our trajectories, once again we bootstrap. Step 5c: Histogram resamples Example directory: examples/AMBER_REMD_with_error_analysis/5-compare_to_raw_canonical Example script: c-histogram_temp_timeseries.sh We can, once again, use multihist.py here. This time to create simple histograms of the canonical resamples. Step 5d: Combine resamples to get distribution with error bars Example directory: examples/AMBER_REMD_with_error_analysis/5-compare_to_raw_canonical Example script: d-combine_resamples_using_avg_dev.sh As with step 4b above. Step 5e: Plot using your favorite plotting program Example directory: examples/AMBER_REMD_with_error_analysis/5-compare_to_raw_canonical I've included an example script (plot.py) using the matplotlib python package. Here's what it looks like:    |