eecr-relay wave shape fit, see the correction to our initial article in JACS. We strongly encourage all the users of the eecr-relay wave shape fit to upgrade, and, unfortunately, refit previously fitted data as the results might change. The other wave shape fits are not affected by the issue.
partially globalparameters using the new
define-indexed-fit command, to pick the best parameters dataset-by-dataset within fit trajectories, but also a parameter space explorer trying all possible permutations of one or more sets of parameters, and the possibility to save the results of a command to a global ruby variable.
There are a lot of other new features, improvements and so on, look for the full list there.
Azide-1.25mm_001.dat, but of course, the purpose of this tutorial is to enable you to work on all of them. The data files contain the time evolution of the absorbance for all wavelengths, in a matrixformat, in which each row correpond to a time point and each column to a wavelength. Start QSoas, and launch the command:
QSoas> load /comments='"'Then, choose the
Azide-1.25mm_001.dat data file. This should bring up a horizontal red line at the bottom of the data display, with X values between about 0 and 2.5. If you zoom on the red line with the mouse wheel, you'll realize it is data. The /comments='"' part is very important since it allows the extraction of the wavelength from the data. We will look at what it means another day. At this stage, you can look at the loaded data using the command:
QSoas> editYou should have a window looking like this:
nan as the X value and what is obviously the wavelength for all the Y values. To tell that QSoas should take this line as the wavelength (which will be the perpendicular coordinate, the coordinate of the other direction of the matrix), first close the edit window and run:
QSoas> set-perp /from-row=0
QSoas> expand /style=red-to-blue /flags=kineticsYour screen should now look like this:
/style=red-to-blue is not strictly necessary, but it gives the red to blue color gradient which makes things easier to look at (and cooler !). The /flags=kinetics is there to attach a label (a flag) to the newly created datasets so we can easily manipulate all of them at the same time. Then it's time to fit, with the following command:
QSoas> mfit-exponential-decay flagged:kineticsThis should bring up a new window. After resizing it, you should have something that looks like this:
Z = 728.15 showing that QSoas has recognized that the currently displayed dataset corresponds to the wavelength 728.15. The equation fitted to the data is: $$y(x) = A_\infty + A_1 \times \exp -(x - x_0)/\tau_1$$
In this case, while the \(A_1\) and \(A_\infty\) parameters clearly depend on the wavelength, the time constant of evolution should be independent of wavelength (the process happens at a certain rate regardless of the wavelength we're analyzing), so that the \(\tau_1\) parameter should be common for all the datasets/wavelengths. Just click on the globalcheckbox at the right of the
tau_1 parameter, make sure it is checked, and hit the Fitbutton... The fit should not take long (less than a minute), and then you end up with the results of the fits: all the parameters. The best way to look at the non global parameters like \(A_1\) and \(A_\infty\) is to use the
Show Parametersitem from the
Parametersmenu. Using it and clicking on
A_inf too should give you a display like this one:
A_inf parameter corresponds to the spectum at infinite time (of azide-bound heme), while the A_1 parameter corresponds to the difference spectrum between the initial (azide-free) and final (azide-bound) states.
Now, the fit is finished, you can save the parameters if you want to reload them in a later fit by using the Parameters/
Savemenu item or export them in a form more suitable for plotting using
Parameters/
Export(although QSoas can also display and the parameters saved using
Save). This concludes this first approach to fitting the data. What you can do is
tau_1 parameter as a function of the azide concentration;QSoas> mfit-exponential-decay /exponentials=2 flagged:kinetics
QSoas> in the code areas above are meant to be typed into the QSoas command line (at the bottom of the window), and started by pressing enter at the end. You must remove the QSoas> bit. The other lines (when applicable) show you the response of QSoas, in the terminal just above the command-line. You may want to play with the QSoas tutorial to learn more about how to interact with QSoas.
kinetic-system fits, and the handling of Marcus-Hush-Chidsey (or Marcus "distribution of states") kinetics for electron transfers.
subfunctionsto run using commands like run-for-each. The second is the possibility to define
variablesto be reused later (like the script arguments) using the new command let. There are a lot of other new features, improvements and so on, look for the full list there.
release branch. I have now decided to publish the currentcode of QSoas in the github repository (in the
public branch). This way, you can follow and use all the good things that were developed since the last release, and also verify whether any bug you have is still present in the currently developed version !eval;split-on-values. The first solution is the one that came naturally to me, and is by far the most general and extensible, but the second one is shorter, and doesn't require external script files.
QSoas> load kcat-vs-ph.dat QSoas> split-on-values pH x /flags=dataAfter these commands, the stacks contains a series of datasets bearing the
data flag, that each contain a single column of data, as can be seen from the beginnings of a show-stack command:
QSoas> k Normal stack: F C Rows Segs Name #0 (*) 1 43 1 'kcat-vs-ph_subset_22.dat' #1 (*) 1 44 1 'kcat-vs-ph_subset_21.dat' #2 (*) 1 43 1 'kcat-vs-ph_subset_20.dat' ...Each of these datasets have a meta-data named
pH whose value is the original x value from kcat-vs-ph.dat. Now, the idea is to run a stats command on the resulting datasets, extracting the average value of x and its standard deviation, together with the value of the meta pH. The most natural and general way to do this is to use run-for-datasets, using the following script file (named process-one.cmds):
stats /meta=pH /output=true /stats=x_average,x_stddevSo the command looks like:
QSoas> run-for-datasets process-one.cmds flagged:dataThis command produces an output file containing, for each flagged dataset, a line containing
x_average, x_stddev, and pH. Then, it is just a matter of loading the output file and shuffling the columns in the right order to get the data in the form asked. Overall, this looks like this:
l kcat-vs-ph.dat split-on-values pH x /flags=data output result.dat /overwrite=true run-for-datasets process-one.cmds flagged:data l result.dat apply-formula tmp=y2;y2=y;y=x;x=tmp dataset-options /yerrors=y2The slight improvement over what is described above is the use of the
output command to write the output to a dedicated file (here result.dat), instead of out.dat and ensuring it is overwritten, so that no data remains from previous runs.
stats command can work with datasets other than the current one, by supplying them to the /buffers= option, so that it is not necessary to use run-for-datasets;l kcat-vs-ph.dat split-on-values pH x /flags=data stats /meta=pH /accumulate=* /stats=x_average,x_stddev /buffers=flagged:data pop apply-formula tmp=y2;y2=y;y=x;x=tmp dataset-options /yerrors=y2
pH value data: the X column is the pH, the Y column the result of the experiment at the given pH (let's say the measure of the catalytic rate of an enzyme). Your task is to take this data and produce a single dataset which contains, for each pH value, the pH, the average of the results at that pH and the standard deviation. The result should be identical to the following file, and should look like that:
stats, and the more natural way in QSoas is to take advantage of split-on-values, which is a very powerful command but somehow hard to master, which is the point of this Quiz.
QSoas> cd QSoas> load 27.oxw QSoas> load 27-blanc.oxw QSoas> S 1 0(after the first command, you have to manually select the directory in which you downloaded the data files). The
S 1 0 command just subtracts the dataset 1 (the first loaded) from the dataset 0 (the last loaded), see more there. blancis the French for
blank... Then, we remove a bit of the beginning and the end of the data, corresponding to one half of the steps at \(E_0\), which we don't exploit much here (they are essentially only used to make sure that the irreversible loss is taken care of properly). This is done using strip-if:
QSoas> strip-if x<30||x>300Then, we can fit ! The fit used is called fit-linear-kinetic-system, which is used to fit kinetic models with only linear reactions (like here) and steps which change the values of the rate constants but do not instantly change the concentrations. The specific command to fit the data is:
QSoas> fit-linear-kinetic-system /species=2 /steps=0,1,2,1,0The
/species=2 indicates that there are two species (A and I). The /steps=0,1,2,1,0 indicates that there are 5 steps, with three different conditions (0 to 2) in order 0,1,2,1,0. This fits needs a bit of setup before getting started. The species are numbered, 1 and 2, and the conditions (potentials) are indicated by #0, #1 and #2 suffixes.
I_1 and I_2 are the currents for the species 1 and 2, so something for 1 (active form) and 0 for 2 (inactive form).
Moreover, the parameters I_2_#0 (and _#1, _#2) should be fixed and not free (since we don't need to adjust a current for the inactive form).k_11 and k_22 correspond to species-specific irreversible loss. It is generally best to leave them fixed to 0.k_12 is the formation of 2 (I) from 1 (A), and k_21 is the formation of A from I. Their values will be determined for the three conditions. The default values should work here.k_loss parameters are the rates of irreversible loss that apply indiscriminately on all species (unlike k_11 and k_22). They are adjusted and ther default values should work too.alpha_1_0 and alpha_2_0 are the initial concentrations of species 1 and 2, so they should be fixed to 1 and 0.xstart_a and (_b, _c, _d and _e) correspond to the starting times for the steps, here, 0, 60, 120, 210 and 270.starting-parameters.params parameters to have all setup the correct way. Then, just hit Fit, enjoy this moment when QSoas works and you don't have to... The screen should now look like this:
21.oxw and 21-blanc.oxw files. Usually, re-loading the best fit parameters from other potentials as starting parameters work really well. Gathering the results of several fits into a real curve of rate constants as a function of potentials is left as an exercise for the reader (or maybe a later post), although you may find these series of posts useful in this context !
fits):
jaggydata.
meta-data. In a previous post, we have already described how one can set meta-data to data that are already loaded, and how one can make use of them. QSoas is already able to figure out some meta-data in the case of electrochemical data, most notably in the case of files acquired by GPES, ECLab or CHI potentiostats. However, only a small number of constructors are supported as of now[1], and there are a number of experimental details that the software is never going to be able to figure out for you, such as the pH, the sample, what you were doing... The new version of QSoas provides a means to permanently store meta-data for experimental data files:
QSoas> record-meta pH 7 file.datThis command uses record-meta to permanently store the information
pH = 7for the file
file.dat. Any time QSoas loads the file again, either today or in one year, the meta-data will contain the value 7 for the field pH. Behind the scenes, QSoas creates a single small file, file.dat.qsm, in which the meta-data are stored (in the form of a JSON dictionnary).
You can set the same meta-data to many files in one go, using wildcards (see load for more information). For instance, to set the pH=7 meta-data to all the .dat files in the current directory, you can use:
QSoas> record-meta pH 7 *.datYou can only set one meta-data for each call to record-meta, but you can use it as many times as you like. Finally, you can use the
/for-which option to load or browse to select only the files which have the meta you need:
QSoas> browse /for-which=$meta.pH<=7This command browses the files in the current directory, showing only the ones that have a pH meta-data which is 7 or below.
[1] I'm always ready to implement the parsing of other file formats that could be useful for you. If you need parsing of special files, please contact me, sending the given files and the meta-data you'd expect to find in those.
$ sign. Using this, we can keep track of the previous values. For instance, to create a new column with the sum of the y values, one can use the following approach:
QSoas> eval $sum=0 QSoas> apply-formula /extra-columns=1 $sum+=y;y2=$sumThe first line initializes the variable to 0, before we start summing, and the code in the second line is run for each dataset row, in order. For the first row, for instance,
$sum is initially 0 (from the eval line); after the execution of the code, it is now the first value of y. After the second row, the second value of y is added, and so on. The image below shows the resulting y2 when used on:
QSoas> generate-dataset -1 1 x
## time current potential 0 0.1 0.5 1 0.2 2 0.3 3 0.2 4 1.2 0.6 5 1.3 ...If you need to have the values everywhere, for instance if you need to split on their values, you could also use a global variable, taking advantage of the fact that missing values are represented by QSoas using "Not A Number" values, which can be detected using the Ruby function
nan?:
QSoas> apply-formula "if y2.nan?; then y2=$value; else $value=y2;end"Note the need of quotes because there are spaces in the ruby code. If the value of
y2 is NaN, that is it is missing, then it is taken from the global variable $value else $value is set the current value of y2. Hence, the values are propagated down:
## time current potential 0 0.1 0.5 1 0.2 0.5 2 0.3 0.5 3 0.2 0.5 4 1.2 0.6 5 1.3 0.6 ...Of course, this doesn't work if the first value of
y2 is missing.
/mode=rename option, does only the renaming part, without saving. You can make full use of meta-data (see also a first post here)for renaming. The full power is unlocked using the /expression= option. For instance, for renaming the last 5 datasets (so numbers 0 to 4) using a scheme based on the value of their pH meta-data, you can use the following code:
QSoas> save-datasets /mode=rename /expression='"dataset-#{$meta.pH}"' 0..4
The double quotes are cumbersome but necessary, since the outer quotes (') prevent the inner ones (") to be removed and the inner quotes are here to indicate to Ruby that we are dealing with text. The bit inside #{...} is interpreted by Ruby as Ruby code; here it is $meta.pH, the value of the "pH" meta-data. Finally the 0..4 specifies the datasets to work with. So theses datasets will change name to become dataset-7 for pH 7, etc...expertmode for fitting
expertmode, turned on by using the
/expert=true option with the fit commands. The expert mode features a command prompt that looks like the normal command prompt, in which it is possible:
parameter space exploration, which consists in running (automatically) a number of fits with different (random) starting parameters to maximize the chances of finding out the best parameters for the fit, avoiding local minima.
Show randomto learn new tricks !).
apply-formula, like automatic code completion;record-meta;set-meta command:
QSoas> set-meta pH 7.5This command sets the meta-data
pH to the value 7.5. Keep in mind that QSoas does not know anything about the meaning of the meta-data[1]. It can keep track of the meta-data you give, and manipulate them, but it will not interpret them for you. You can set several meta-data by repeating calls to set-meta, and you can display the meta-data attached to a dataset using the command show. Here is an example:
QSoas> generate-buffer 0 10 QSoas> set-meta pH 7.5 QSoas> set-meta sample "My sample" QSoas> show 0 Dataset generated.dat: 2 cols, 1000 rows, 1 segments, #0 Flags: Meta-data: pH = 7.5 sample = My sampleNote here the use of quotes around
My sample since there is a space inside the value.
1 is run, a line with the information about the largest peak of the current dataset is written to the output file. It is possible to automatically add meta-data to those lines by using the /meta= option of the output command. Just listing the names of the meta-data will add them to each line of the output file.
As a full example, we'll see how one can take advantage of meta-data to determine the position of the peak of the function \(x^2 \exp (-a\,x)\) depends on \(a\). For that, we first create a script that generates the function for a certain value of \(a\), sets the meta-data a to the corresponding value, and find the peak. Let's call this file do-one.cmds (all the script files can be found in the GitHub repository):
generate-buffer 0 20 x**2*exp(-x*${1})
set-meta a ${1}
1
This script takes a single argument, the value of \(a\), generates the appropriate dataset, sets the meta-data a and writes the data about the largest (and only in this case) peak to the output file. Let's now run this script with 1 as an argument:
QSoas> @ do-one.cmds 1This command generates a file
out.dat containing the following data:
## buffer what x y index width left_width right_width area generated.dat max 2.002002002 0.541340590883 100 3.4034034034 1.24124124124 2.162162162161.99999908761This gives various information about the peak found: the name of the dataset it was found in, whether it's a maximum or minimum, the x and y positions of the peak, the index in the file, the widths of the peak and its area. We are interested here mainly in the
x position.
Then, we just run this script for several values of \(a\) using run-for-each, and in particular the option /range-type=lin that makes it interpret values like 0.5..5:80 as 80 values evenly spread between 0.5 and 5. The script is called
run-all.cmds:
output peaks.dat /overwrite=true /meta=a run-for-each do-one.cmds /range-type=lin 0.5..5:80 V all /style=red-to-blueThe first line sets up the output to the output file
peaks.dat. The option /meta=a makes sure the meta a is added to each line of the output file, and /overwrite=true make sure the file is overwritten just before the first data is written to it, in order to avoid accumulating the results of different runs of the script. The last line just displays all the curves with a color gradient. It looks like this:
@ run-all.cmds) creates a new file peaks.dat, whose first line looks like this:
## buffer what x y index width left_width right_width area aThe column
x (the 3rd) contains the position of the peaks, and the column a (the 10th) contains the meta a (this column wasn't present in the output we described above, because we had not used yet the output /meta=a command). Therefore, to load the peak position as a function of a, one has just to run:
QSoas> load peaks.dat /columns=10,3This looks like this:
[1] this is not exactly true. For instance, some commands like unwrap interpret the sr meta-data as a voltammetric scan rate if it is present. But this is the exception.
load-one.cmds:
load ${1}
apply-formula y2=y**2 /extra-columns=1
flag /flags=processed
When this script is run with the name of a spectrum file as argument, it loads it (replaces ${1} by the first argument, the file name), adds a column y2 containing the square of the y column, and flag it with the processed flag. This is not absolutely necessary, but it makes it much easier to refer to all the spectra when they are processed. Then to process all the spectra, one just has to run the following commands:
run-for-each load-one.cmds Spectrum-1.dat Spectrum-2.dat Spectrum-3.dat average flagged:processed apply-formula y2=(y2-y**2)**0.5 dataset-options /yerrors=y2The
run-for-each command runs the load-one.cmds script for all the spectra (one could also have used Spectra-*.dat to not have to give all the file names). Then, the average averages the values of the columns over all the datasets. To be clear, it finds all the values that have the same X (or very close X values) and average them, column by column. The result of this command is therefore a dataset with the average of the original \(y\) data as y column and the average of the original \(y^2\) data as y2 column. So now, the only thing left to do is to use the above equation, which is done by the apply-formula code. The last command, dataset-options, is not absolutely necessary but it signals to QSoas that the standard error of the y column should be found in the y2 column. This is now available as script method-one.cmds in the git repository.
apply-formula, the value $stats.y_stddev corresponds to the standard deviation of the whole y column... Loading the spectra yields just a series of x,y datasets. We can contract them into a single dataset with one x column and several y columns:
load Spectrum-*.dat /flags=spectra contract flagged:spectraAfter these commands, the current dataset contains data in the form of:
lambda1 a1_1 a1_2 a1_3 lambda2 a2_1 a2_2 a2_3 ...in which the
ai_1 come from the first file, ai_2 the second and so on. We need to use transpose to transform that dataset into:
0 a1_1 a2_1 ... 1 a1_2 a2_2 ... 2 a1_3 a2_3 ...In this dataset, values of the absorbance for the same wavelength for each dataset is now stored in columns. The next step is just to use
expand to obtain a series of datasets with the same x column and a single y column (each corresponding to a different wavelength in the original data). The game is now to replace these datasets with something that looks like:
0 a_average 1 a_stddevFor that, one takes advantage of the
$stats.y_average and $stats.y_stddev values in apply-formula, together with the i special variable that represents the index of the point:
apply-formula "if i == 0; then y=$stats.y_average; end; if i == 1; then y=$stats.y_stddev; end" strip-if i>1Then, all that is left is to apply this to all the datasets created by
expand, which can be just made using run-for-datasets, and then, we reverse the splitting by using contract and transpose ! In summary, this looks like this. We need two files. The first, process-one.cmds contains the following code:
apply-formula "if i == 0; then y=$stats.y_average; end; if i == 1; then y=$stats.y_stddev; end" strip-if i>1 flag /flags=processedThe main file,
method-two.cmds looks like this:
load Spectrum-*.dat /flags=spectra contract flagged:spectra transpose expand /flags=tmp run-for-datasets process-one.cmds flagged:tmp contract flagged:processed transpose dataset-options /yerrors=y2Note some of the code above can be greatly simplified using new features present in the upcoming 3.0 version, but that is the topic for another post.
fit-arb fit:
QSoas> l michaelis.dat QSoas> fit-arb vm/(1+km/x)After running the fit, the window should look like this:
vm) and \(K_m\) (km). But, for publication, one would want to generate "smooth" curve going through the lines... Saving the curve from "Data.../Save all" doesn't help, since the data has as many points as the original data and looks very "jaggy" (like on the screenshot above)... So one needs a curve with more data points.
Maybe the most natural solution is simply to use generate-buffer together with apply-formula using the formula and the values of km and vm obtained from the fit, like:
QSoas> generate-buffer 0 20 QSoas> apply-formula y=3.51742/(1+3.69767/x)By default,
generate-buffer generate 1000 evenly spaced x values, but you can change their number using the /samples option. The two above commands can be combined to just one call to generate-buffer:
QSoas> generate-buffer 0 20 3.51742/(1+3.69767/x)This works, but it is quite cumbersome and it is not going to work well for complex formulas or the results of differential equations or kinetic systems... This is why to each
fit- command corresponds a sim- command that computes the result of the fit using a "saved parameters" file (here, michaelis.params, but you can also save it yourself) and buffers as "models" for X values:
QSoas> generate-buffer 0 20 QSoas> sim-arb vm/(1+km/x) michaelis.params 0This strategy works with every single fit ! As an added benefit, you even get the fit parameters as meta-data, which are displayed by the
show command:
QSoas> show 0 Dataset generated_fit_arb.dat: 2 cols, 1000 rows, 1 segments, #0 Flags: Meta-data: commands = sim-arb vm/(1+km/x) michaelis.params 0 fit = arb (formula: vm/(1+km/x)) km = 3.69767 vm = 3.5174They also get saved as comments if you save the data. Important note: the
sim-arb command will be available only in the 3.0 release, although you can already enjoy it if you use the github version.
average command:
QSoas> load Spectrum*.dat /flags=spectra QSoas> average flagged:spectraHowever,
average does not provide means to make standard deviations, it just takes the average of all but the X column. I wanted to add this feature, but I realized there are already at least two distinct ways to do that...Spectrum-1.dat, Spectrum-2.dat and Spectrum-3.dat). There are at least two ways:
average and on apply-formula, and which requires that you remember how to compute standard deviations.contract for instance) and relies on the fact that you can use statistics in apply-formula (and in particular you can use y_stddev to refer to the standard deviation of \(y\)), but which does not require you to know exactly how to compute standard deviations.Average.dat. The figure below shows a zoom on the data superimposed to the average (bonus points to find how to display this light red area that corresponds to the standard deviation !).
To begin this tutorial, first download the files from the github repository (direct links: data, parameter file and ruby script). Start QSoas, go to the directory where you saved the files, load the data file, and remove spikes in the data using the following commands:
QSoas> cd QSoas> l Km-CODH-IV.dat QSoas> R
QSoas> fit-arb im/(1+km/s) /with=s:1,expThis simply launches the fit interface to fit the exact equations above. The
im/(1+km/s) is simply the translation of the Michaelis-Menten equation above, and the /with=s:1,exp specifies that s is the result of the sum of 1 exponential like for the definition of above. Then, load the Km-CODH-IV.params parameter files (using the "Parameters.../Load from file" action at the bottom, or the Ctrl+L keyboard shortcut). Your window should now look like this:
Ctrl+F).
QSoas> fit-arb im/(1+km/s)+io /with=s:1,expNotice the
+io bit that corresponds to the addition of an offset current. Load again the base parameters, run the fit again... Your fit window show now look like:
im is \(i_m\), the maximal current, around 120 nA (which matches the magnitude of the original current).io is the offset current, around -3nA.km is the \(K_m\), expressed in the same units as s_1, the first "injected" value of s (we used 50 because the injection is 50 µM CO). That means the value of \(K_m\) determined by this fit is around 9 nM ! (but see below).s_t_1 is the time of the injection of CO (it was in the parameter files you loaded).s_tau_1 is the time constant of departure of CO, noted \(\tau\) in the equations above.
transport.rb Ruby file that contains the definition of the itrpt function, and re-launch the fit window using:
QSoas> ruby-run transport.rb QSoas> fit-arb itrprt(s,km,nFAm,nFAmu)+io /with=s:1,expLoad again the parameter file... but this time you'll have to play a bit more with the starting parameters for QSoas to find the right values when you fit. Here are some tips:
io, no parameter should be negative;
1.0. They are also freely available for download from the downloads page.ruby...ruby end block. Probably the most useful possibility is to define elaborated functions directly from within QSoas, or, preferable, from within a script; this is an alternative to defining a function in a completely separated Ruby-only file using ruby-run. For instance, you can define a function for plain Michaelis-Menten kinetics with a file containing:ruby def my_func(x, vm, km) return vm/(1 + km/x) end ruby end
my_func with three parameters, , (vm) and (km), with the formula: QSoas> eval my_func(1.0,1.0,1.0) => 0.5 QSoas> eval my_func(1e4,1.0,1.0) => 0.999900009999
x = 1.0, vm = 1.0 and km = 1.0. For , the result is (here 0.5). For , the result is almost . You can use the newly defined my_func in any place you would use any ruby code, such as in the optional argument to generate-buffer, or for arbitrary fits:QSoas> generate-buffer 0 10 my_func(x,3.0,0.6) QSoas> fit-arb my_func(x,vm,km)
my_func, just run the ruby code again with a new definition, such as:ruby def my_func(x, vm, km) return vm/(1 + km/x**2) end ruby endThe previous version is just erased, and all new uses of
my_func will refer to your new definition.qsoas-goodies github repository for more goodies !