Re: [Ifeffit] R-factor uncertainty
Hi Lisa,
At 07:23 AM 1/5/2007, you wrote:
I have a question about the R factor: how can I decide if the difference between the R factors of 2 fits is statistically significant, i.e, how can I calculate the uncertainty which has to be associated to the R factor?
As I understand it, you can't. The R-factor is not a proper statistical measure, as it doesn't incorporate any measure of data quality. That's the great weakness of this measure of quality-of-fit. It is also its strength, as estimating the uncertainty in EXAFS data is notoriously problematic. The complementary statistic is reduced chi-square. It does incorporate a measure of data quality. By default, ifeffit uses noise from high in the FT to estimate this. That's a reasonable idea, but can be problematic. It has been shown (by Matt and/or Shelly, as I recall), that there may in some cases be signal in the part of the FT ifeffit is using to estimate noise. There are also cases where the noise may not be "white," that is, the noise high in the FT may be a poor estimate of the noise low in the FT. Ifeffit does allow you to specify a value of the measurement uncertainty instead, so if you think you have a way of doing this, go ahead. What does all this mean in practice? It means, in my opinion, that the actual =value= of the reduced chi-square statistic is usually meaningless, unless you have a good way of coming up with the measurement uncertainty (for example, your sample may be so dilute that errors are dominated by counting statistics). But reduced chi-square is a great statistic for comparing two fits to a given set of data, particularly if the k-range, k-weighting, and k-window are the same for the two fits. For example, you can apply statistical tests of significance, if you'd like. The R-factor then provides the reality check that the fit is "good" at all. The R-factor isn't doing anything other than what you can see by looking at a graph, but is a nice shorthand for tables showing the results of many fits and similar applications. If there's a big R-factor (say, 0.20), the question of statistical significance isn't necessary to tell you that you haven't got a conclusive positive result: maybe the R-factor is big because the fitting model is lousy, or maybe it's big because the data quality is lousy, but either way the fit shouldn't be trusted. I'd also add that your eye tells you considerably more than the R-factor, because you can tell the character of the mismatch. Is it in the high part of the FT, low, or evenly throughout? Is the miss primarily in amplitude, or phase? I often find I choose a fit with an R-factor of 0.03 over one with 0.01, if, for example, the 0.03 reproduces qualitatively all the features in the data but has small errors in the amplitude of the peaks, while the 0.01 fits the first part of the spectrum perfectly but misses some peak altogether. Hope that helps... --Scott Calvin Sarah Lawrence College
Hello Artemis users, I found an interesting feature in Artemis when generating a list of atoms from the attached cif file worth sharing on the list. In the attached cif file for Pd2Si the position of the Si2 atoms are at x,y,z = 0.3333, 0.6667, 0 with a space group of p -6 2 m. Using this cif file and choosing either Pd1 or Pd2 as the core atom results in errors in the atom list. For example for the Pd2 site there are 4 Si neighbors with the 1st and 3rd and the 2nd and 4th on top of each other. If you edit the atomic positions so that the Si2 atoms are correct to an additional digit x,y,z = 0.33333, 0.66667,0 then the atom list is fixed with only two Pd2-Si neighbors at each of the unique positions. The interesting feature is that the coordinates on the atoms page need to be correct to the 5 decimal place. Cheers, Shelly
On Saturday 06 January 2007 10:53, Kelly, Shelly D. wrote:
Hello Artemis users,
I found an interesting feature in Artemis when generating a list of atoms from the attached cif file worth sharing on the list. In the attached cif file for Pd2Si the position of the Si2 atoms are at x,y,z = 0.3333, 0.6667, 0 with a space group of p -6 2 m. Using this cif file and choosing either Pd1 or Pd2 as the core atom results in errors in the atom list. For example for the Pd2 site there are 4 Si neighbors with the 1st and 3rd and the 2nd and 4th on top of each other. If you edit the atomic positions so that the Si2 atoms are correct to an additional digit x,y,z = 0.33333, 0.66667,0 then the atom list is fixed with only two Pd2-Si neighbors at each of the unique positions. The interesting feature is that the coordinates on the atoms page need to be correct to the 5 decimal place.
This has been discussed many, many times on the mailing list and is explained in question 3 at http://cars9.uchicago.edu/iffwiki/FAQ/RunningFeff The version of Feff that ships with Ifeffit (and, I believe, all others) uses 5 significant figures in the path finder. Since Atoms main purpose is to interact with Feff, it interprets the 5 digit thing quite strictly. B -- Bruce Ravel ---------------------------------------------- bravel@anl.gov Molecular Environmental Science Group, Building 203, Room E-165 MRCAT, Sector 10, Advanced Photon Source, Building 433, Room B007 Argonne National Laboratory phone and voice mail: (1) 630 252 5033 Argonne IL 60439, USA fax: (1) 630 252 9793 My homepage: http://cars9.uchicago.edu/~ravel EXAFS software: http://cars9.uchicago.edu/~ravel/software/
Hi all, I have a little problem using ifeffit with program code taken from athena. I try to make some linear combination calculations. Using Athena it works but using the code with ifeffit, I get some errors. Everything works fine untill the following lines: unguess erase @group l___cf set l___cf.eee = rh.energy+0 set l___cf_npts = npts(l___cf.energy) set l___cf.data = splint(l___cf.eee, rh.flat, l___cf.energy) ** cannot set scalar value to array ** l___cf.data = splint(l___cf.eee,rh.flat,l___cf.energy) The problem here seems to be the array "l___cf.energy". It is empty and I never found some declarations for it. At the beginning of my program I also load all macros athena is loading. Someone have an idea what is wrong with that code? Thank you for help, best regards Bernd P.S. If it helps, here is the complete code for ifeffit of my perl program: # load macros ifeffit(" load(/root/.ifeffit/macros.iff)"); #get standards ifeffit(" read_data(file = /home/promo/Cu/pics/test/Cu.txt.spek, group=\"cu\", label = 'energy_KeV xmu')"); ifeffit(" set cu.energy = cu.energy_KeV * 1000"); ifeffit(" set ___n = npts(cu.energy_KeV)"); ifeffit(" pre_edge(\"cu.energy+0\", cu.xmu, e0=8979.1, pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430)"); ifeffit(" spline(\"cu.energy+0\", cu.xmu, e0=8979.1, rbkg=1.0, kmin=0.5, kmax=12, kweight=1, dk=2, kwindow=kaiser-bessel, \ pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430, interp=quad)"); ifeffit(" set cu.preline = pre_offset + pre_slope * (cu.energy + 0)"); ifeffit(" set cu.postline = norm_c0 + norm_c1 * (cu.energy + 0.000000) + norm_c2 * (cu.energy + 0.000000)**2"); ifeffit(" step cu.energy e0 cu.theta"); ifeffit(" set flat_c0 = norm_c0 - pre_offset"); ifeffit(" set flat_c1 = norm_c1 - pre_slope"); ifeffit(" set flat_c2 = norm_c2"); ifeffit(" set cu.line = (flat_c0 + flat_c1 * (cu.energy+0) + flat_c2 * (cu.energy+0)**2)"); ifeffit(" set cu.flat = (cu.pre + (edge_step - cu.line) * cu.theta) / edge_step"); ifeffit(" set cu.fbkg = (cu.bkg-cu.preline+(edge_step-cu.line)*cu.theta)/edge_step"); ifeffit(" read_data(file = /home/promo/Cu/pics/test/CuO.txt.spek, group=\"cuo\", label = 'energy_KeV xmu')"); ifeffit(" set cuo.energy = cuo.energy_KeV * 1000"); ifeffit(" set ___n = npts(cuo.energy_KeV)"); ifeffit(" pre_edge(\"cuo.energy+0\", cuo.xmu, e0=8991.6, pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430)"); ifeffit(" spline(\"cuo.energy+0\", cuo.xmu, e0=8991.6, rbkg=1.0, kmin=0.5, kmax=12, kweight=1, dk=2, kwindow=kaiser-bessel, \ pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430, interp=quad)"); ifeffit(" set cuo.preline = pre_offset + pre_slope * (cuo.energy + 0)"); ifeffit(" set cuo.postline = norm_c0 + norm_c1 * (cuo.energy + 0.000000) + norm_c2 * (cuo.energy + 0.000000)**2"); ifeffit(" step cuo.energy e0 cuo.theta"); ifeffit(" set flat_c0 = norm_c0 - pre_offset"); ifeffit(" set flat_c1 = norm_c1 - pre_slope"); ifeffit(" set flat_c2 = norm_c2"); ifeffit(" set cuo.line = (flat_c0 + flat_c1 * (cuo.energy+0) + flat_c2 * (cuo.energy+0)**2)"); ifeffit(" set cuo.flat = (cuo.pre + (edge_step - cuo.line) * cuo.theta) / edge_step"); ifeffit(" set cuo.fbkg = (cuo.bkg-cuo.preline+(edge_step-cuo.line)*cuo.theta)/edge_step"); ifeffit(" read_data(file = /home/promo/Cu/pics/test/Cu2O.txt.spek, group=\"cu2o\", label = 'energy_KeV xmu')"); ifeffit(" set cu2o.energy = cu2o.energy_KeV * 1000"); ifeffit(" set ___n = npts(cu2o.energy_KeV)"); ifeffit(" pre_edge(\"cu2o.energy+0\", cu2o.xmu, e0=8980.6, pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430)"); ifeffit(" spline(\"cu2o.energy+0\", cu2o.xmu, e0=8980.6, rbkg=1.0, kmin=0.5, kmax=12, kweight=1, dk=2, kwindow=kaiser-bessel, \ pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430, interp=quad)"); ifeffit(" set cu2o.preline = pre_offset + pre_slope * (cu2o.energy + 0)"); ifeffit(" set cu2o.postline = norm_c0 + norm_c1 * (cu2o.energy + 0.000000) + norm_c2 * (cu2o.energy + 0.000000)**2"); ifeffit(" step cu2o.energy e0 cu2o.theta"); ifeffit(" set flat_c0 = norm_c0 - pre_offset"); ifeffit(" set flat_c1 = norm_c1 - pre_slope"); ifeffit(" set flat_c2 = norm_c2"); ifeffit(" set cu2o.line = (flat_c0 + flat_c1 * (cu2o.energy+0) + flat_c2 * (cu2o.energy+0)**2)"); ifeffit(" set cu2o.flat = (cu2o.pre + (edge_step - cu2o.line) * cu2o.theta) / edge_step"); ifeffit(" set cu2o.fbkg = (cu2o.bkg-cu2o.preline+(edge_step-cu2o.line)*cu2o.theta)/edge_step"); #get data ifeffit(" read_data(file = /tmp/spektrum$i.dat, group=\"rh\", label = 'energy_KeV xmu')"); ifeffit(" set rh.energy = rh.energy_KeV * 1000"); ifeffit(" set ___n = npts(rh.energy_KeV)"); ifeffit(" pre_edge(\"rh.energy+0\", rh.xmu, pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430)"); ifeffit(" spline(\"rh.energy+0\", rh.xmu, rbkg=1.0, kmin=0.5, kmax=12, kweight=1, dk=2, kwindow=kaiser-bessel, \ pre1=-150, pre2=-30, norm_order=3, norm1=150, norm2=430, interp=quad)"); $estart = get_scalar("e0"); ifeffit(" set rh.preline = pre_offset + pre_slope * (rh.energy + 0)"); ifeffit(" set rh.postline = norm_c0 + norm_c1 * (rh.energy + 0.000000) + norm_c2 * (rh.energy + 0.000000)**2"); ifeffit(" step rh.energy e0 rh.theta"); ifeffit(" set flat_c0 = norm_c0 - pre_offset"); ifeffit(" set flat_c1 = norm_c1 - pre_slope"); ifeffit(" set flat_c2 = norm_c2"); ifeffit(" set rh.line = (flat_c0 + flat_c1 * (rh.energy+0) + flat_c2 * (rh.energy+0)**2)"); ifeffit(" set rh.flat = (rh.pre + (edge_step - rh.line) * rh.theta) / edge_step"); ifeffit(" set rh.fbkg = (rh.bkg-rh.preline+(edge_step-rh.line)*rh.theta)/edge_step"); ifeffit(" set &status = 0"); ifeffit(" erase delta_e1 delta_w1 delta_ww1"); ifeffit(" erase delta_e2 delta_w2 delta_ww2"); ifeffit(" erase delta_e3 delta_w3 delta_ww3"); ifeffit(" erase delta_e4 delta_w4 delta_ww4"); ifeffit(" erase delta_yint delta_slope"); # performing linear combination fit in norm(E) ifeffit(" unguess"); ifeffit(" erase \@group l___cf"); ifeffit(" set l___cf.eee = rh.energy+0"); ifeffit(" set l___cf_npts = npts(l___cf.energy)"); ifeffit(" set l___cf.data = splint(l___cf.eee, rh.flat, l___cf.energy)"); ifeffit(" guess e1 = 0"); ifeffit(" guess ww1 = 0.333"); ifeffit(" def w1 = max(0,min(ww1,1))"); ifeffit(" def l___cf.1 = abs(w1)*splint(cu.energy+e1+0.00, cu.flat, l___cf.energy)"); ifeffit(" def e2 = e1"); ifeffit(" guess ww2 = 0.333"); ifeffit(" def w2 = max(0,min(ww2,1))"); ifeffit(" def l___cf.2 = abs(w2)*splint(cuo.energy+e2+0.00, cuo.flat, l___cf.energy)"); ifeffit(" def w3 = max(0, 1 - (w1+w2))"); ifeffit(" def e3 = e1"); ifeffit(" def l___cf.3 = w3*splint(cu20.energy+e3+0.00, cu2o.flat, l___cf.energy)"); ifeffit(" def l___cf.mix = l___cf.1 + l___cf.2 + l___cf.3"); ifeffit(" def l___cf.resid = l___cf.mix - l___cf.data"); $xmin = $estart - 20; $xmax = $estart + 30; ifeffit(" minimize(l___cf.resid, x=l___cf.energy, $xmin, $xmax)"); ...
On Sunday 07 January 2007 20:31, Bernd Griesebock wrote:
I have a little problem using ifeffit with program code taken from athena. I try to make some linear combination calculations. Using Athena it works but using the code with ifeffit, I get some errors. Everything works fine untill the following lines:
unguess erase @group l___cf set l___cf.eee = rh.energy+0 set l___cf_npts = npts(l___cf.energy) set l___cf.data = splint(l___cf.eee, rh.flat, l___cf.energy)
** cannot set scalar value to array ** l___cf.data = splint(l___cf.eee,rh.flat,l___cf.energy)
The problem here seems to be the array "l___cf.energy". It is empty and I never found some declarations for it. At the beginning of my program I also load all macros athena is loading. Someone have an idea what is wrong with that code?
For some reason, when I wrote the code that generates those lines of Ifeffit commands, I decided not to use Ifefift's slice and nofx commands to construct the energy array. The thing that is missing from that sequence of commands is the bit where the energy array is set to the contents of the eee array between emin and emax. This is done because the interpolations of the data and the minimization are then done on the grid of the fit. Rather than doing that chore using Ifeffit commands, I chose to manipulate those arrays directly in perl. I wrote that code almost 3 years ago, so I don't remember my reasoning. The lines I suggest below should do what you need done, assuming that you replace emin and emax with the correct values of absolute energy of the range over which you are doing the fit. set l___cf.eee = rh.energy+0 set xmin = nofx(l___cf.eee, emin) set xmax = nofx(l___cf.eee, emax) set l___cf.energy = slice(l___cf.eee, xmin, xmax) set l___cf_npts = npts(l___cf.energy) set l___cf.data = splint(l___cf.eee, rh.flat, l___cf.energy) B -- Bruce Ravel ---------------------------------------------- bravel@anl.gov Molecular Environmental Science Group, Building 203, Room E-165 MRCAT, Sector 10, Advance Photon Source, Building 433, Room B007 Argonne National Laboratory phone and voice mail: (1) 630 252 5033 Argonne IL 60439, USA fax: (1) 630 252 9793 My homepage: http://cars9.uchicago.edu/~ravel EXAFS software: http://cars9.uchicago.edu/~ravel/software/
Thank you very much Bruce, that help me a lot. I finally fixed it by using a slightly modulated part of the source code of Athena. At this point I just want to report a little bug I found using the linear combination function of Athena. Since I didn't use the plot window I didn't realized it before. initial position: opened in Athena for example 3 standards and 2 groups of data sets. Highlighted group1, switching to linear combination fitting (lcf). If one switch now to group2 (still in lcf-window) the results of the fit wont be what one expects cause the xmin and xmax for minimization function is not recalculated. There can be a big difference in results if one uses for example 2 groups of data sets, where e0 values are not really close to each other. I hope that maybe helps a bit. Best regards, Bernd Bruce Ravel schrieb:
On Sunday 07 January 2007 20:31, Bernd Griesebock wrote:
I have a little problem using ifeffit with program code taken from athena. I try to make some linear combination calculations. Using Athena it works but using the code with ifeffit, I get some errors. Everything works fine untill the following lines:
unguess erase @group l___cf set l___cf.eee = rh.energy+0 set l___cf_npts = npts(l___cf.energy) set l___cf.data = splint(l___cf.eee, rh.flat, l___cf.energy)
** cannot set scalar value to array ** l___cf.data = splint(l___cf.eee,rh.flat,l___cf.energy)
The problem here seems to be the array "l___cf.energy". It is empty and I never found some declarations for it. At the beginning of my program I also load all macros athena is loading. Someone have an idea what is wrong with that code?
For some reason, when I wrote the code that generates those lines of Ifeffit commands, I decided not to use Ifefift's slice and nofx commands to construct the energy array. The thing that is missing from that sequence of commands is the bit where the energy array is set to the contents of the eee array between emin and emax. This is done because the interpolations of the data and the minimization are then done on the grid of the fit. Rather than doing that chore using Ifeffit commands, I chose to manipulate those arrays directly in perl. I wrote that code almost 3 years ago, so I don't remember my reasoning.
The lines I suggest below should do what you need done, assuming that you replace emin and emax with the correct values of absolute energy of the range over which you are doing the fit.
set l___cf.eee = rh.energy+0
set xmin = nofx(l___cf.eee, emin) set xmax = nofx(l___cf.eee, emax) set l___cf.energy = slice(l___cf.eee, xmin, xmax)
set l___cf_npts = npts(l___cf.energy) set l___cf.data = splint(l___cf.eee, rh.flat, l___cf.energy)
B
participants (4)
-
Bernd Griesebock -
Bruce Ravel -
Kelly, Shelly D. -
Scott Calvin