diff --git a/qiskit_experiments/curve_analysis/__init__.py b/qiskit_experiments/curve_analysis/__init__.py index d08f11a989..3115d395e3 100644 --- a/qiskit_experiments/curve_analysis/__init__.py +++ b/qiskit_experiments/curve_analysis/__init__.py @@ -467,6 +467,7 @@ def _create_analysis_results(self, fit_data, quality, **metadata): BaseCurveAnalysis CurveAnalysis + MultiGroupCurveAnalysis Data Classes ============ @@ -548,6 +549,7 @@ def _create_analysis_results(self, fit_data, quality, **metadata): """ from .base_curve_analysis import BaseCurveAnalysis from .curve_analysis import CurveAnalysis +from .grouped_curve_analysis import MultiGroupCurveAnalysis from .curve_data import ( CurveData, CurveFitResult, diff --git a/qiskit_experiments/curve_analysis/base_curve_analysis.py b/qiskit_experiments/curve_analysis/base_curve_analysis.py index bd7b12aa65..cda50376bf 100644 --- a/qiskit_experiments/curve_analysis/base_curve_analysis.py +++ b/qiskit_experiments/curve_analysis/base_curve_analysis.py @@ -86,8 +86,6 @@ class BaseCurveAnalysis(BaseAnalysis, ABC): This method to creates analysis results for important fit parameters that might be defined by analysis options ``result_parameters``. - In addition, another entry for all fit parameters is created when - the analysis option ``return_fit_parameters`` is ``True``. .. rubric:: _create_curve_data @@ -391,6 +389,7 @@ def _run_curve_fit( self, curve_data: CurveData, models: List[lmfit.Model], + init_guesses: Dict[str, float], ) -> CurveFitResult: """Perform curve fitting on given data collection and fit models. @@ -398,6 +397,8 @@ def _run_curve_fit( curve_data: Formatted data to fit. models: A list of LMFIT models that are used to build a cost function for the LMFIT minimizer. + init_guesses: Dictionary of fit parameter initial guesses keyed on the + fit parameter name. Returns: The best fitting outcome with minimum reduced chi-squared value. @@ -419,10 +420,11 @@ def _run_curve_fit( default_fit_opt = FitOptions( parameters=unite_parameter_names, - default_p0=self.options.p0, + default_p0=init_guesses, default_bounds=self.options.bounds, **self.options.lmfit_options, ) + # Bind fixed parameters if not empty if self.options.fixed_parameters: fixed_parameters = { diff --git a/qiskit_experiments/curve_analysis/curve_analysis.py b/qiskit_experiments/curve_analysis/curve_analysis.py index 9db69780d3..9b277c13fc 100644 --- a/qiskit_experiments/curve_analysis/curve_analysis.py +++ b/qiskit_experiments/curve_analysis/curve_analysis.py @@ -193,6 +193,7 @@ def _run_analysis( fit_data = self._run_curve_fit( curve_data=formatted_data, models=self._models, + init_guesses=self.options.p0, ) if fit_data.success: @@ -201,15 +202,15 @@ def _run_analysis( quality = "bad" if self.options.return_fit_parameters: - # Store fit status entry regardless of success. + # Store fit status overview entry regardless of success. # This is sometime useful when debugging the fitting code. - fit_parameters = AnalysisResultData( + overview = AnalysisResultData( name=PARAMS_ENTRY_PREFIX + self.__class__.__name__, value=fit_data, quality=quality, extra=self.options.extra, ) - analysis_results.append(fit_parameters) + analysis_results.append(overview) # Create figure and result data if fit_data.success: diff --git a/qiskit_experiments/curve_analysis/grouped_curve_analysis.py b/qiskit_experiments/curve_analysis/grouped_curve_analysis.py new file mode 100644 index 0000000000..9169137d5e --- /dev/null +++ b/qiskit_experiments/curve_analysis/grouped_curve_analysis.py @@ -0,0 +1,324 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2021. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +""" +Analysis class for multi-group curve fitting. +""" +# pylint: disable=invalid-name + +from typing import Dict, List, Tuple, Optional, Any +from collections import defaultdict + +import lmfit +import numpy as np +from uncertainties import unumpy as unp, UFloat + +from qiskit_experiments.framework import ExperimentData, AnalysisResultData, Options +from qiskit_experiments.exceptions import AnalysisError + +from .base_curve_analysis import BaseCurveAnalysis, PARAMS_ENTRY_PREFIX +from .curve_data import CurveFitResult +from .utils import analysis_result_to_repr, eval_with_uncertainties + + +class MultiGroupCurveAnalysis(BaseCurveAnalysis): + r"""Curve analysis with multiple independent groups. + + The :class:`.MultiGroupCurveAnalysis` provides a fitting framework for the case + + .. math:: + + \Theta_{i, \mbox{opt}} = \arg\min_{\Theta_{i, \mbox{fit}}} (F(X_i, \Theta_i)-Y_i)^2, + + where :math:`F` is a common fit model for multiple independent dataset :math:`(X_i, Y_i)`. + The fit generates multiple fit outcomes :math:`\Theta_{i, \mbox{opt}}` for each dataset. + + Each experiment circuit must have metadata of data sorting keys for model mapping and + group mapping. For example, if the experiment has two experiment conditions "A" and "B" + while fitting the outcomes with three models "X", "Y", "Z", then the metadata must consist + of the condition and model, e.g. ``{"condition": "A", "model": "X"}``. + + In this example, "condition" should be defined in the ``group_data_sort_key`` + of the :class:`.MultiGroupCurveAnalysis`, and "model" should be defined in + the ``data_sort_key`` for each LMFIT model. + + .. code-block:: python + + import lmfit + from qiskit_experiments.curve_analysis import MultiGroupCurveAnalysis + + analysis = MultiGroupCurveAnalysis( + groups=["groupA", "groupB"], + group_data_sort_key=[{"condition": "A"}, {"condition": "B"}], + models=[ + lmfit.Model(fit_func_x, name="curve_X", data_sort_key={"model": "X"}), + lmfit.Model(fit_func_y, name="curve_Y", data_sort_key={"model": "Y"}), + lmfit.Model(fit_func_z, name="curve_Z", data_sort_key={"model": "Z"}), + ] + ) + + In above setting, each fit curve will appear in the output figure with unique name + in the form of ``{group name}_{model name}``, e.g. ``groupA_curve_X``. + + A subclass can implement the following extra method. + + .. rubric:: _create_composite_analysis_results + + This method computes new quantities by taking all fit outcomes and + store the new values in the analysis results. + By default, no extra quantity is computed. + + """ + + def __init__( + self, + groups: List[str], + group_data_sort_key: List[Dict[str, Any]], + models: Optional[List[lmfit.Model]] = None, + ): + """Create new analysis. + + Args: + groups: List of group names. + group_data_sort_key: List of sort key dictionary for each group. + models: List of LMFIT ``Model`` class to define fitting functions and + parameters. Provided models are shared among groups. + + Raises: + AnalysisError: When groups and group_data_sort_key don't match. + """ + super().__init__() + + if len(groups) != len(group_data_sort_key): + raise AnalysisError( + "Number of groups and data sort key doesn't match. " + "Data sort key must be provided for each group." + ) + + self._groups = groups + self._group_data_sort_key = group_data_sort_key + self._models = models + + @classmethod + def _default_options(cls) -> Options: + """Default analysis options. + + Analysis Options: + p0 (Union[List[Dict[str, float]], Dict[str, float]]): A dictionary of initial guesses + for the fit parameters keyed on the fit parameter names. + This can be provided for either each fit group or entire groups. + To provide initial guesses for each group, a list of dictionary must be + set here. The list ordering corresponds to the groups. + """ + return super()._default_options() + + @property + def parameters(self) -> List[str]: + """Return parameters of this curve analysis.""" + unite_params = [] + for model in self._models: + for name in model.param_names: + if name not in unite_params and name not in self.options.fixed_parameters: + unite_params.append(name) + return unite_params + + @property + def models(self) -> List[lmfit.Model]: + """Return fit models.""" + return self._models + + # pylint: disable=unused-argument + def _create_composite_analysis_results( + self, + fit_dataset: List[CurveFitResult], + quality: str, + **metadata, + ) -> List[AnalysisResultData]: + """Create new analysis data from combination of all fit outcomes. + + Args: + fit_dataset: Fit outcomes. The ordering of data corresponds to the + groups defined in the ``self._groups``. + quality: Overall quality of fittings. + + Returns: + List of analysis results. + """ + return [] + + def _run_analysis( + self, + experiment_data: ExperimentData, + ) -> Tuple[List[AnalysisResultData], List["matplotlib.figure.Figure"]]: + + self._initialize(experiment_data) + analysis_results = [] + + # Extract experiment data for this group + # TODO use dataframe + group_data = defaultdict(list) + for datum in experiment_data.data(): + for group, dkey in zip(self._groups, self._group_data_sort_key): + if all(datum["metadata"].get(k, None) == v for k, v in dkey.items()): + group_data[group].append(datum) + break + + fit_qualities = [] + fit_dataset = [] + for group in self._groups: + metadata = self.options.extra.copy() + metadata["group"] = group + + # Run data processing + processed_data = self._run_data_processing( + raw_data=group_data[group], + models=self._models, + ) + + if self.options.plot and self.options.plot_raw_data: + for model in self._models: + sub_data = processed_data.get_subset_of(model._name) + self.drawer.draw_raw_data( + x_data=sub_data.x, + y_data=sub_data.y, + name=f"{model._name}_{group}", + ) + + # Format data + formatted_data = self._format_data(processed_data) + if self.options.plot: + for model in self._models: + sub_data = formatted_data.get_subset_of(model._name) + self.drawer.draw_formatted_data( + x_data=sub_data.x, + y_data=sub_data.y, + y_err_data=sub_data.y_err, + name=f"{model._name}_{group}", + ) + + # Run fitting + if isinstance(self.options.p0, dict): + initial_guesses = self.options.p0 + else: + initial_guesses = self.options.p0[self._groups.index(group)] + + fit_data = self._run_curve_fit( + curve_data=formatted_data, + models=self._models, + init_guesses=initial_guesses, + ) + + if fit_data.success: + quality = self._evaluate_quality(fit_data) + else: + quality = "bad" + + if self.options.return_fit_parameters: + overview = AnalysisResultData( + name=f"{PARAMS_ENTRY_PREFIX}{self.__class__.__name__}_{group}", + value=fit_data, + quality=quality, + extra=metadata, + ) + analysis_results.append(overview) + + if fit_data.success: + # Add extra analysis results + analysis_results.extend( + self._create_analysis_results( + fit_data=fit_data, quality=quality, **metadata.copy() + ) + ) + + # Draw fit result + if self.options.plot: + interp_x = np.linspace( + np.min(formatted_data.x), np.max(formatted_data.x), num=100 + ) + for model in self._models: + y_data_with_uncertainty = eval_with_uncertainties( + x=interp_x, + model=model, + params=fit_data.ufloat_params, + ) + y_mean = unp.nominal_values(y_data_with_uncertainty) + # Draw fit line + self.drawer.draw_fit_line( + x_data=interp_x, + y_data=y_mean, + name=f"{model._name}_{group}", + ) + if fit_data.covar is not None: + # Draw confidence intervals with different n_sigma + sigmas = unp.std_devs(y_data_with_uncertainty) + if np.isfinite(sigmas).all(): + for n_sigma, alpha in self.drawer.options.plot_sigma: + self.drawer.draw_confidence_interval( + x_data=interp_x, + y_ub=y_mean + n_sigma * sigmas, + y_lb=y_mean - n_sigma * sigmas, + name=f"{model._name}_{group}", + alpha=alpha, + ) + + # Add raw data points + analysis_results.extend( + self._create_curve_data(curve_data=formatted_data, models=self._models) + ) + + fit_qualities.append(quality) + fit_dataset.append(fit_data) + + quality = "good" if all(q == "good" for q in fit_qualities) else "bad" + + # Create analysis results by combining all fit data + analysis_results.extend( + self._create_composite_analysis_results( + fit_dataset=fit_dataset, quality=quality, **self.options.extra.copy() + ) + ) + + if self.options.plot: + # Write fitting report + report = "" + for res in analysis_results: + if isinstance(res.value, (float, UFloat)): + report += f"{analysis_result_to_repr(res)}\n" + chisqs = [] + for group, fit_data in zip(self._groups, fit_dataset): + chisqs.append(r"reduced-$\chi^2$ = " + f"{fit_data.reduced_chisq: .4g} ({group})") + report += "\n".join(chisqs) + self.drawer.draw_fit_report(description=report) + + # Finalize plot + if self.options.plot: + self.drawer.format_canvas() + return analysis_results, [self.drawer.figure] + + return analysis_results, [] + + def __getstate__(self): + state = self.__dict__.copy() + # Convert models into JSON str. + # This object includes local function and cannot be pickled. + source = [m.dumps() for m in state["_models"]] + state["_models"] = source + return state + + def __setstate__(self, state): + model_objs = [] + for source in state.pop("_models"): + tmp_mod = lmfit.Model(func=None) + mod = tmp_mod.loads(s=source) + model_objs.append(mod) + self.__dict__.update(state) + self._models = model_objs diff --git a/qiskit_experiments/library/characterization/analysis/cr_hamiltonian_analysis.py b/qiskit_experiments/library/characterization/analysis/cr_hamiltonian_analysis.py index d85efb2031..2114a3957d 100644 --- a/qiskit_experiments/library/characterization/analysis/cr_hamiltonian_analysis.py +++ b/qiskit_experiments/library/characterization/analysis/cr_hamiltonian_analysis.py @@ -12,8 +12,6 @@ """Cross resonance Hamiltonian tomography experiment analysis.""" -from collections import defaultdict -from itertools import product from typing import List, Union import lmfit @@ -25,8 +23,7 @@ from qiskit_experiments.framework import AnalysisResultData -# pylint: disable=line-too-long -class CrossResonanceHamiltonianAnalysis(curve.CurveAnalysis): +class CrossResonanceHamiltonianAnalysis(curve.MultiGroupCurveAnalysis): r"""A class to analyze cross resonance Hamiltonian tomography experiment. # section: fit_model @@ -36,35 +33,35 @@ class CrossResonanceHamiltonianAnalysis(curve.CurveAnalysis): .. math:: \begin{align} - F_{x, c}(t) &= \frac{1}{\Omega_c^2} \left( - - p_{z, c} p_{x, c} + p_{z, c} p_{x, c} \cos(\Omega_c t') + - \Omega_c p_{y, c} \sin(\Omega_c t') \right) + b \tag{1} \\ - F_{y, c}(t) &= \frac{1}{\Omega_c^2} \left( - p_{z, c} p_{y, c} - p_{z, c} p_{y, c} \cos(\Omega_c t') - - \Omega_c p_{x, c} \sin(\Omega_c t') \right) + b \tag{2} \\ - F_{z, c}(t) &= \frac{1}{\Omega_c^2} \left( - p_{z, c}^2 + (p_{x, c}^2 + p_{y, c}^2) \cos(\Omega_c t') \right) + b \tag{3} + F_x(t) &= \frac{1}{\Omega^2} \left( + - p_z p_x + p_z p_x \cos(\Omega t') + + \Omega p_y \sin(\Omega t') \right) + b \tag{1} \\ + F_y(t) &= \frac{1}{\Omega^2} \left( + p_z p_y - p_z p_y \cos(\Omega t') - + \Omega p_x \sin(\Omega t') \right) + b \tag{2} \\ + F_z(t) &= \frac{1}{\Omega^2} \left( + p_z^2 + (p_x^2 + p_y^2) \cos(\Omega t') \right) + b \tag{3} \end{align} where :math:`t' = t + t_{\rm offset}` with :math:`t` is pulse duration to scan and :math:`t_{\rm offset}` is an extra fit parameter that may represent the edge effect. - The :math:`\Omega_c = \sqrt{p_{x, c}^2+p_{y, c}^2+p_{z, c}^2}` and - :math:`p_{x, c}, p_{y, c}, p_{z, c}, b` are also fit parameters. - The subscript :math:`c` represents the state of control qubit :math:`c \in \{0, 1\}`. - The fit functions :math:`F_{x, c}, F_{y, c}, F_{z, c}` approximate the Pauli expectation - values :math:`\langle \sigma_{x, c} (t) \rangle, \langle \sigma_{y, c} (t) \rangle, - \langle \sigma_{z, c} (t) \rangle` of the target qubit, respectively. + The :math:`\Omega = \sqrt{p_x^2+p_y^2+p_z^2}`, and :math:`p_x, p_y, p_z, b` are fit parameters. + The fit functions :math:`F_x, F_y, F_z` approximate the Pauli expectation + values :math:`\langle \sigma_x (t) \rangle, \langle \sigma_y (t) \rangle, + \langle \sigma_z (t) \rangle` of the target qubit, respectively. + This fitting is performed for two experiment results for different control qubit state + in :\math:`c \in [|0\rangle, |1\rangle]`. Based on the fit result, cross resonance Hamiltonian coefficients can be written as .. math:: - ZX &= \frac{p_{x, 0} - p_{x, 1}}{2} \\ - ZY &= \frac{p_{y, 0} - p_{y, 1}}{2} \\ - ZZ &= \frac{p_{z, 0} - p_{z, 1}}{2} \\ - IX &= \frac{p_{x, 0} + p_{x, 1}}{2} \\ - IY &= \frac{p_{y, 0} + p_{y, 1}}{2} \\ - IZ &= \frac{p_{z, 0} + p_{z, 1}}{2} + ZX &= \frac{p_{x, |0\rangle} - p_{x, |1\rangle}}{2} \\ + ZY &= \frac{p_{y, |0\rangle} - p_{y, |1\rangle}}{2} \\ + ZZ &= \frac{p_{z, |0\rangle} - p_{z, |1\rangle}}{2} \\ + IX &= \frac{p_{x, |0\rangle} + p_{x, |1\rangle}}{2} \\ + IY &= \frac{p_{y, |0\rangle} + p_{y, |1\rangle}}{2} \\ + IZ &= \frac{p_{z, |0\rangle} + p_{z, |1\rangle}}{2} In this analysis, the initial guess is generated by the following equations. @@ -89,33 +86,18 @@ class CrossResonanceHamiltonianAnalysis(curve.CurveAnalysis): .parametric_pulses.GaussianSquare` pulse envelope. bounds: [0, None] - defpar p_{x, 0}: - desc: Fit parameter of oscillations when control qubit state is 0. + defpar p_x: + desc: Fit parameter of oscillations. init_guess: See fit model section. bounds: None - defpar p_{y, 0}: - desc: Fit parameter of oscillations when control qubit state is 0. + defpar p_y: + desc: Fit parameter of oscillations. init_guess: See fit model section. bounds: None - defpar p_{z, 0}: - desc: Fit parameter of oscillations when control qubit state is 0. - init_guess: See fit model section. - bounds: None - - defpar p_{x, 1}: - desc: Fit parameter of oscillations when control qubit state is 1. - init_guess: See fit model section. - bounds: None - - defpar p_{y, 1}: - desc: Fit parameter of oscillations when control qubit state is 1. - init_guess: See fit model section. - bounds: None - - defpar p_{z, 1}: - desc: Fit parameter of oscillations when control qubit state is 1. + defpar p_z: + desc: Fit parameter of oscillations. init_guess: See fit model section. bounds: None @@ -139,22 +121,23 @@ def __init__(self): } models = [] - for state in (0, 1): - for axis, temp_eq in eqr_temps.items(): - eq = temp_eq - for op in ("px", "py", "pz"): - eq = eq.replace(op, f"{op}{state}") - eq = eq.replace("W", f"sqrt(px{state}**2 + py{state}**2 + pz{state}**2)") - eq = eq.replace("X", "(x + t_off)") - models.append( - lmfit.models.ExpressionModel( - expr=eq, - name=f"{axis}|c={state}", - data_sort_key={"control_state": state, "meas_basis": axis}, - ) + for axis, temp_eq in eqr_temps.items(): + eq = temp_eq + eq = eq.replace("W", "sqrt(px**2 + py**2 + pz**2)") + eq = eq.replace("X", "(x + t_off)") + models.append( + lmfit.models.ExpressionModel( + expr=eq, + name=f"{axis}", + data_sort_key={"meas_basis": axis}, ) + ) - super().__init__(models=models) + super().__init__( + groups=["ctrl0", "ctrl1"], + group_data_sort_key=[{"control_state": 0}, {"control_state": 1}], + models=models, + ) @classmethod def _default_options(cls): @@ -174,12 +157,12 @@ def _default_options(cls): fit_report_rpos=(0.28, -0.10), ylim=(-1, 1), plot_options={ - "x|c=0": {"color": "blue", "symbol": "o", "canvas": 0}, - "y|c=0": {"color": "blue", "symbol": "o", "canvas": 1}, - "z|c=0": {"color": "blue", "symbol": "o", "canvas": 2}, - "x|c=1": {"color": "red", "symbol": "^", "canvas": 0}, - "y|c=1": {"color": "red", "symbol": "^", "canvas": 1}, - "z|c=1": {"color": "red", "symbol": "^", "canvas": 2}, + "x_ctrl0": {"color": "blue", "symbol": "o", "canvas": 0}, + "y_ctrl0": {"color": "blue", "symbol": "o", "canvas": 1}, + "z_ctrl0": {"color": "blue", "symbol": "o", "canvas": 2}, + "x_ctrl1": {"color": "red", "symbol": "^", "canvas": 0}, + "y_ctrl1": {"color": "red", "symbol": "^", "canvas": 1}, + "z_ctrl1": {"color": "red", "symbol": "^", "canvas": 2}, }, ) default_options.data_processor = dp.DataProcessor( @@ -203,86 +186,63 @@ def _generate_fit_guesses( Returns: List of fit options that are passed to the fitter function. """ + fit_options = [] + user_opt.bounds.set_if_empty(t_off=(0, np.inf), b=(-1, 1)) user_opt.p0.set_if_empty(b=1e-9) - guesses = defaultdict(list) - for control in (0, 1): - x_data = curve_data.get_subset_of(f"x|c={control}") - y_data = curve_data.get_subset_of(f"y|c={control}") - z_data = curve_data.get_subset_of(f"z|c={control}") - - omega_xyz = [] - for data in (x_data, y_data, z_data): - ymin, ymax = np.percentile(data.y, [10, 90]) - if ymax - ymin < 0.2: - # oscillation amplitude might be almost zero, - # then exclude from average because of lower SNR - continue - fft_freq = curve.guess.frequency(data.x, data.y) - omega_xyz.append(fft_freq) - if omega_xyz: - omega = 2 * np.pi * np.average(omega_xyz) - else: - omega = 1e-3 - - zmin, zmax = np.percentile(z_data.y, [10, 90]) - theta = np.arccos(np.sqrt((zmax - zmin) / 2)) - - # The FFT might be up to 1/2 bin off - df = 2 * np.pi / ((z_data.x[1] - z_data.x[0]) * len(z_data.x)) - for omega_shifted in [omega, omega - df / 2, omega + df / 2]: - for phi in np.linspace(-np.pi, np.pi, 5): - px = omega_shifted * np.cos(theta) * np.cos(phi) - py = omega_shifted * np.cos(theta) * np.sin(phi) - pz = omega_shifted * np.sin(theta) - guesses[control].append( - { - f"px{control}": px, - f"py{control}": py, - f"pz{control}": pz, - } - ) - if omega < df: - # empirical guess for low frequency case - guesses[control].append( - { - f"px{control}": omega, - f"py{control}": omega, - f"pz{control}": 0, - } + x_data = curve_data.get_subset_of("x") + y_data = curve_data.get_subset_of("y") + z_data = curve_data.get_subset_of("z") + + omega_xyz = [] + for data in (x_data, y_data, z_data): + ymin, ymax = np.percentile(data.y, [10, 90]) + if ymax - ymin < 0.2: + # oscillation amplitude might be almost zero, + # then exclude from average because of lower SNR + continue + fft_freq = curve.guess.frequency(data.x, data.y) + omega_xyz.append(fft_freq) + if omega_xyz: + omega = 2 * np.pi * np.average(omega_xyz) + else: + omega = 1e-3 + + zmin, zmax = np.percentile(z_data.y, [10, 90]) + theta = np.arccos(np.sqrt((zmax - zmin) / 2)) + + # The FFT might be up to 1/2 bin off + df = 2 * np.pi / ((z_data.x[1] - z_data.x[0]) * len(z_data.x)) + for omega_shifted in [omega, omega - df / 2, omega + df / 2]: + for phi in np.linspace(-np.pi, np.pi, 5): + new_opt = user_opt.copy() + new_opt.p0.set_if_empty( + px=omega_shifted * np.cos(theta) * np.cos(phi), + py=omega_shifted * np.cos(theta) * np.sin(phi), + pz=omega_shifted * np.sin(theta), ) - - fit_options = [] - # combine all guesses in Cartesian product - for p0s, p1s in product(guesses[0], guesses[1]): + fit_options.append(new_opt) + if omega < df: + # empirical guess for low frequency case new_opt = user_opt.copy() - new_opt.p0.set_if_empty(**p0s, **p1s) + new_opt.p0.set_if_empty(px=omega, py=omega, pz=0) fit_options.append(new_opt) return fit_options - def _create_analysis_results( + def _create_composite_analysis_results( self, - fit_data: curve.CurveFitResult, + fit_dataset: List[curve.CurveFitResult], quality: str, **metadata, ) -> List[AnalysisResultData]: - """Create analysis results for important fit parameters. - - Args: - fit_data: Fit outcome. - quality: Quality of fit outcome. - - Returns: - List of analysis result data. - """ - outcomes = super()._create_analysis_results(fit_data, quality, **metadata) + outcomes = [] for control in ("z", "i"): for target in ("x", "y", "z"): - p0_val = fit_data.ufloat_params[f"p{target}0"] - p1_val = fit_data.ufloat_params[f"p{target}1"] + p0_val = fit_dataset[0].ufloat_params[f"p{target}"] + p1_val = fit_dataset[1].ufloat_params[f"p{target}"] if control == "z": coef_val = 0.5 * (p0_val - p1_val) / (2 * np.pi) @@ -293,7 +253,6 @@ def _create_analysis_results( AnalysisResultData( name=f"omega_{control}{target}", value=coef_val, - chisq=fit_data.reduced_chisq, quality=quality, extra={ "unit": "Hz", diff --git a/releasenotes/notes/curve-analysis-02a702a81e014adf.yaml b/releasenotes/notes/curve-analysis-02a702a81e014adf.yaml index c79cb86f2a..e930b08e13 100644 --- a/releasenotes/notes/curve-analysis-02a702a81e014adf.yaml +++ b/releasenotes/notes/curve-analysis-02a702a81e014adf.yaml @@ -11,6 +11,14 @@ features: The default fit method is ``"least_squares"``. Analysis class author can flexibly define new analysis instance with LMFIT ``Model`` objects. See LMFIT documentation for user guide. + - | + New curve analysis baseclass :class:`.MultiGroupCurveAnalysis` has been added. + This curve analysis variant offers a framework to fit multiple indenepdent datasets + with the consistent fit model. For example, if you define an experiment scanning + a parameter with different conditions, e.g. with different control qubit states + in some two-qubit gate experiment, the multi-group curve analysis can implement + its analysis with simpler code comparing with writing a batch experment. + See class documentation for more details. - | New options have been added to to the :class:`.CurveAnaysis` curve drawer. diff --git a/test/curve_analysis/test_baseclass.py b/test/curve_analysis/test_baseclass.py index 605cdfb42c..15239542e5 100644 --- a/test/curve_analysis/test_baseclass.py +++ b/test/curve_analysis/test_baseclass.py @@ -21,7 +21,7 @@ from lmfit.models import ExpressionModel from qiskit.qobj.utils import MeasLevel -from qiskit_experiments.curve_analysis import CurveAnalysis, fit_function +from qiskit_experiments.curve_analysis import CurveAnalysis, MultiGroupCurveAnalysis, fit_function from qiskit_experiments.curve_analysis.curve_data import ( SeriesDef, CurveFitResult, @@ -457,6 +457,64 @@ def test_get_init_params(self): y_reproduced = analysis.models[0].eval(x=x, **overview.init_params) np.testing.assert_array_almost_equal(y_ref, y_reproduced) + def test_multi_group_curve_analysis(self): + """Integration test for multi group curve analysis.""" + + analysis = MultiGroupCurveAnalysis( + groups=["group_A", "group_B"], + group_data_sort_key=[{"cond": 0}, {"cond": 1}], + models=[ + ExpressionModel( + expr="amp * cos(2 * pi * freq * x) + b", + data_sort_key={"type": "cos"}, + ), + ExpressionModel( + expr="amp * sin(2 * pi * freq * x) + b", + data_sort_key={"type": "sin"}, + ), + ], + ) + analysis.set_options( + p0=[{"amp": 0.3, "freq": 2.1, "b": 0.5}, {"amp": 0.5, "freq": 3.2, "b": 0.5}], + result_parameters=["amp"], + plot=False, + ) + + amp1 = 0.2 + amp2 = 0.4 + b1 = 0.5 + b2 = 0.5 + freq1 = 2.1 + freq2 = 3.2 + + x = np.linspace(0, 1, 100) + y1a = amp1 * np.cos(2 * np.pi * freq1 * x) + b1 + y2a = amp1 * np.sin(2 * np.pi * freq1 * x) + b1 + y1b = amp2 * np.cos(2 * np.pi * freq2 * x) + b2 + y2b = amp2 * np.sin(2 * np.pi * freq2 * x) + b2 + + test_data1a = self.single_sampler(x, y1a, type="cos", cond=0) + test_data2a = self.single_sampler(x, y2a, type="sin", cond=0) + test_data1b = self.single_sampler(x, y1b, type="cos", cond=1) + test_data2b = self.single_sampler(x, y2b, type="sin", cond=1) + + expdata = ExperimentData(experiment=FakeExperiment()) + expdata.add_data(test_data1a.data()) + expdata.add_data(test_data2a.data()) + expdata.add_data(test_data1b.data()) + expdata.add_data(test_data2b.data()) + expdata.metadata["meas_level"] = MeasLevel.CLASSIFIED + + result = analysis.run(expdata).block_for_results() + amps = result.analysis_results("amp") + + # two entries are generated for group A and group B + self.assertEqual(len(amps), 2) + self.assertEqual(amps[0].extra["group"], "group_A") + self.assertEqual(amps[1].extra["group"], "group_B") + self.assertAlmostEqual(amps[0].value.n, 0.2, delta=0.1) + self.assertAlmostEqual(amps[1].value.n, 0.4, delta=0.1) + class TestFitOptions(QiskitExperimentsTestCase): """Unittest for fit option object."""