Code reference

Module level aliases

For user convenience, the following objects are available at the module level.

class nanite.Indentation

alias of nanite.indent.Indentation

class nanite.IndentationGroup

alias of nanite.group.IndentationGroup

class nanite.IndentationRater

alias of nanite.rate.IndentationRater

class nanite.QMap

alias of nanite.qmap.QMap

nanite.load_group()

alias of nanite.group.load_group

Force-indentation data

class nanite.indent.Indentation(data, metadata, diskcache=None)[source]

Additional functionalities for afmformats.AFMForceDistance

apply_preprocessing(preprocessing=None, options=None, ret_details=False)[source]

Perform curve preprocessing steps

Parameters:
  • preprocessing (list) – A list of preprocessing method identifiers that are stored in the nanite.preproc.PREPROCESSORS list. If set to None, self.preprocessing will be used.

  • options (dict of dict) – Dictionary of keyword arguments for each preprocessing step (if applicable)

  • ret_details – Return preprocessing details dictionary

compute_emodulus_mindelta(callback=None)[source]

Elastic modulus in dependency of maximum indentation

The fitting interval is varied such that the maximum indentation depth ranges from the lowest tip position to the estimated contact point. For each interval, the current model is fitted and the elastic modulus is extracted.

Parameters:

callback (callable) – A method that is called with the emoduli and indentations as the computation proceeds every five steps.

Returns:

emoduli, indentations – The fitted elastic moduli at the corresponding maximal indentation depths.

Return type:

1d ndarrays

Notes

The information about emodulus and mindelta is also stored in self.fit_properties with the keys “optimal_fit_E_array” and “optimal_fit_delta_array”, if self.fit_model is called with the argument search_optimal_fit set to True.

estimate_contact_point_index(method='deviation_from_baseline')[source]

Estimate the contact point index

See the poc submodule for more information.

estimate_optimal_mindelta()[source]

Estimate the optimal indentation depth

This is a convenience function that wraps around compute_emodulus_mindelta and IndentationFitter.compute_opt_mindelta.

fit_model(**kwargs)[source]

Fit the approach-retract data to a model function

Parameters:
  • model_key (str) – A key referring to a model in nanite.model.models_available

  • params_initial (instance of lmfit.Parameters or dict) – Parameters for fitting. If not given, default parameters are used.

  • range_x (tuple of 2) – The range for fitting, see range_type below.

  • range_type (str) –

    One of:

    • absolute:

      Set the absolute fitting range in values given by the x_axis.

    • relative cp:

      In some cases it is desired to be able to fit a model only up until a certain indentation depth (tip position) measured from the contact point. Since the contact point is a fit parameter as well, this requires a two-pass fitting.

  • preprocessing (list of str) – Preprocessing

  • preprocessing_options (list of str) – Preprocessing

  • segment (str) – Segment index (e.g. 0 for approach)

  • weight_cp (float) – Weight the contact point region which shows artifacts that are difficult to model with e.g. Hertz.

  • gcf_k (float) – Geometrical correction factor \(k\) for non-single-contact data. The measured indentation is multiplied by this factor to correct for experimental geometries during fitting, e.g. gcf_k=0.5 for parallel-place compression.

  • optimal_fit_edelta (bool) – Search for the optimal fit by varying the maximal indentation depth and determining a plateau in the resulting Young’s modulus (fitting parameter “E”).

get_ancillary_parameters(model_key=None)[source]

Compute ancillary parameters for the current model

get_initial_fit_parameters(model_key=None, common_ancillaries=True, model_ancillaries=True)[source]

Return the initial fit parameters

If there are not initial fit parameters set in self.fit_properties, then they are computed.

Parameters:
  • model_key (str) – Optionally set a model key. This will override the “model_key” key in self.fit_properties.

  • common_ancillaries (bool) – Guess global ancillaries such as the contact point.

  • model_ancillaries (bool) – Guess model-related ancillaries

Notes

global_ancillaries and model_ancillaries only have an effect if self.fit_properties[“params_initial”] is set.

get_rating_parameters()[source]

Return current rating parameters

rate_quality(regressor='Extra Trees', training_set='zef18', names=None, lda=None)[source]

Compute the quality of the obtained curve

Uses heuristic approaches to rate a curve.

Parameters:
  • regressor (str) – The regressor name used for rating.

  • training_set (str) – A label for a training set shipped with nanite or a path to a training set.

  • names (list of str) – Only use these features for rating

  • lda (bool) – Perform linear discriminant analysis

Returns:

rating – A value between 0 and 10 where 0 is the lowest rating. If no fit has been performed, a rating of -1 is returned.

Return type:

float

Notes

The rating is cached based on the fitting hash (see IndentationFitter._hash).

property data
property fit_properties

Fitting results, see Indentation.fit_model())

preprocessing

Default preprocessing steps, see Indentation.apply_preprocessing().

preprocessing_options

Preprocessing options

Groups

class nanite.group.IndentationGroup(path=None, meta_override=None, callback=None)[source]

Group of Indentation

Parameters:
  • path (str or pathlib.Path or None) – The path to the data file. The data format is determined and the file is loaded using index.

  • meta_override (dict) – if specified, contains key-value pairs of metadata that should be used when loading the files (see afmformats.meta.META_FIELDS)

  • callback (callable or None) – A method that accepts a float between 0 and 1 to externally track the process of loading the data.

append(afmdata)[source]

Append a new instance of AFMData

This subclassed method makes sure that “spring constant” is set if “tip position” has to be computed in the future.

Parameters:

afmdata (afmformats.afm_data.AFMData) – AFM data

nanite.group.load_group(path, callback=None, meta_override=None)[source]

Load indentation data from disk

Parameters:
  • path (path-like) – Path to experimental data

  • callback (callable) – function for tracking progress; must accept a float in [0, 1] as an argument.

  • meta_override (dict) – if specified, contains key-value pairs of metadata that should be used when loading the files (see afmformats.meta.META_FIELDS)

Returns:

group – Indentation group with force-distance data

Return type:

nanite.IndetationGroup

Loading data

nanite.read.get_data_paths(path)[source]

Return list of data paths with force-distance data

DEPRECATED

nanite.read.get_data_paths_enum(path, skip_errors=False)[source]

Return a list with paths and their internal enumeration

Parameters:
  • path (str or pathlib.Path or list of str or list of pathlib.Path) – path to data files or directory containing data files; if directories are given, they are searched recursively

  • skip_errors (bool) – skip paths that raise errors

Returns:

path_enum – each entry in the list is a list of [pathlib.Path, int], enumerating all curves in each file

Return type:

list of lists

nanite.read.get_load_data_modality_kwargs()[source]

Return imaging modality kwargs for afmformats.load_data

Uses DEFAULT_MODALITY.

Returns:

kwargs – keyword arguments for afmformats.load_data()

Return type:

dict

nanite.read.load_data(path, callback=None, meta_override=None)[source]

Load data and return list of afmformats.AFMForceDistance

This is essentially a wrapper around afmformats.formats.find_data() and afmformats.formats.load_data() that returns force-distance datasets.

Parameters:
  • path (str or pathlib.Path or list of str or list of pathlib.Path) – path to data files or directory containing data files; if directories are given, they are searched recursively

  • callback (callable) – function for progress tracking; must accept a float in [0, 1] as an argument.

  • meta_override (dict) – if specified, contains key-value pairs of metadata that are used when loading the files (see afmformats.meta.META_FIELDS)

nanite.read.DEFAULT_MODALITY = 'force-distance'

The default imaging modality when loading AFM data. Set this to None to also be able to load e.g. creep-compliance data. See issue https://github.com/AFM-analysis/nanite/issues/11 for more information. Note that especially the export of rating containers may not work with any imaging modality other than force-distance.

Preprocessing

exception nanite.preproc.CannotSplitWarning[source]
class nanite.preproc.IndentationPreprocessor[source]
apply(**kwargs)
autosort(**kwargs)
available(**kwargs)
check_order(**kwargs)
get_func(**kwargs)
get_name(**kwargs)
get_steps_required(**kwargs)
nanite.preproc.apply(apret, identifiers=None, options=None, ret_details=False, preproc_names=None)[source]

Perform force-distance preprocessing steps

Parameters:
  • apret (nanite.Indentation) – The afm data to preprocess

  • identifiers (list) – A list of preprocessing identifiers that will be applied (in the order given).

  • options (dict of dict) – Preprocessing options for each identifier

  • ret_details – Return preprocessing details dictionary

  • preproc_names (list) – Deprecated - use identifiers instead

nanite.preproc.autosort(identifiers)[source]

Automatically sort preprocessing identifiers

This takes into account steps_required and steps_optional.

nanite.preproc.available()[source]

Return list of available preprocessor identifiers

nanite.preproc.check_order(identifiers)[source]

Check preprocessing steps for correct order

nanite.preproc.find_turning_point(tip_position, force, contact_point_index)[source]

Compute the turning point for a force-indentation curve

The turning point is defined as the point that is farthest away from the contact point in the direction of indentation. This implementation normalizes tip position according to their minimum and maximum values. This is necessary, because they live on different orders of magnitudes/units.

nanite.preproc.get_func(identifier)[source]

Return preprocessor function for identifier

nanite.preproc.get_name(identifier)[source]

Return preprocessor name for identifier

nanite.preproc.get_steps_required(identifier)[source]

Return requirement identifiers for identifier

nanite.preproc.preproc_compute_tip_position(apret)[source]

Perform tip-sample separation

Populate the “tip position” column by adding the force normalized by the spring constant to the cantilever height (“height (measured)”).

This computation correctly reproduces the column “Vertical Tip Position” as it is exported by the JPK analysis software with the checked option “Use Unsmoothed Height”.

nanite.preproc.preproc_correct_force_offset(apret)[source]

Correct the force offset with an average baseline value

nanite.preproc.preproc_correct_force_slope(apret, region='baseline', strategy='shift', ret_details=False)[source]

Subtract a linear slope from selected parts of the force curve

Slope correction is a debatable topic in AFM analysis. Many voices are of the opinion that this falsifies data and should not be done. But when data are scarce, slope correction can help to extract valuable information.

Slope correction uses the baseline (the part of the curve before the POC) to determine the slope that should be subtracted from the data. Note that for tilted data, you should use a contact point estimate based on a line and a polynomial.

There are three regions for curve correction:

  • baseline: Only the data leading up to the contact point are modified. This has the advantage that indentation data are not modified, but the baseline can still be used (as weight) in the fitting scheme.

  • appraoch: The slope is subtracted from the entire approach curve. This makes sense when you know that the slope affects baseline and indentation part of your dataset.

  • all: The entire dataset is corrected. This makes sense if you would like to extract information about e.g. viscosity from the area between the approach and retract curves of your dataset.

In addition, there are two strategies available:

  • The “drift” approach assumes that the there is a global drift in the dataset (caused e.g. by a thermal drift). Since these things are usually of temporal nature, the correction in done over the time axis. Use this if the beggining of the approach part and the end of the retract part of your dataset “stick out” in opposite directions when plotting force over tip position.

  • The “shift” approach assumes that there is a global shift that is a function of the distance between sample and cantilever. Use this if the beginning of the approach part and the end of the retract part align when plotting force over tip position.

nanite.preproc.preproc_correct_split_approach_retract(apret)[source]

Split the approach and retract curves (farthest point method)

Approach and retract curves are defined by the microscope. When the direction of piezo movement is flipped, the force at the sample tip is still increasing. This can be either due to a time lag in the AFM system or due to a residual force acting on the sample due to the bent cantilever.

To repair this time lag, we append parts of the retract curve to the approach curve, such that the curves are split at the minimum height.

nanite.preproc.preproc_correct_tip_offset(apret, method='deviation_from_baseline', ret_details=False)[source]

Estimate the point of contact

An estimate of the contact point is subtracted from the tip position.

nanite.preproc.preproc_smooth_height(apret)[source]

Make height data monotonic

For the columns “height (measured)”, “height (piezo), and “tip position”, this method ensures that the approach and retract segments are monotonic.

nanite.preproc.preprocessing_step(identifier, name, steps_required=None, steps_optional=None, options=None)[source]

Decorator for Indentation preprocessors

The name and identifier are stored as a property of the wrapped function.

Parameters:
  • identifier (str) – identifier of the preprocessor (e.g. “correct_tip_offset”)

  • name (str) – human-readble name of the preprocessor (e.g. “Estimate contact point”)

  • steps_required (list of str) – list of preprocessing steps that must be added before this step

  • steps_optional (list of str) – unlike steps_required, these steps do not have to be set, but if they are set, they should come before this step

  • options (list of dict) – if the preprocessor accepts optional keyword arguments, this list yields valid values or dtypes

nanite.preproc.PREPROCESSORS = [<function preproc_compute_tip_position>, <function preproc_correct_force_offset>, <function preproc_correct_force_slope>, <function preproc_correct_tip_offset>, <function preproc_correct_split_approach_retract>, <function preproc_smooth_height>]

Available preprocessors

Contact point estimation

Methods for estimating the point of contact (POC)

nanite.poc.compute_poc(force, method='deviation_from_baseline', ret_details=False)[source]

Compute the contact point from force data

Parameters:
  • force (1d ndarray) – Force data

  • method (str) – Name of the method for computing the POC (see POC_METHODS)

  • ret_details (bool) – Whether or not to return a dictionary with details alongside the POC estimate.

Notes

If the POC method returns np.nan, then the center of the force data is returned (to allow fitting algorithms to proceed).

nanite.poc.compute_preproc_clip_approach(force)[source]

Clip the approach part (discard the retract part)

This POC preprocessing method may be applied before applying the POC estimation method.

nanite.poc.poc(identifier, name, preprocessing)[source]

Decorator for point of contact (POC) methods

The name and identifier are stored as a property of the wrapped function.

Parameters:
  • identifier (str) – identifier of the POC method (e.g. “baseline_deviation”)

  • name (str) – human-readble name of the POC method (e.g. “Deviation from baseline”)

  • preprocessing (list of str) – list of preprocessing methods that should be applied; may contain [“clip_approach”].

nanite.poc.poc_deviation_from_baseline(force, ret_details=False)[source]

Deviation from baseline

  1. Obtain the baseline (initial 10% of the gradient curve)

  2. Compute average and maximum deviation of the baseline

  3. The CP is the index of the curve where it exceeds twice of the maximum deviation

nanite.poc.poc_fit_constant_line(force, ret_details=False)[source]

Piecewise fit with constant and line

Fit a piecewise function (constant+linear) to the baseline and indentation part:

\[F = \text{max}(d, m\delta + d)\]

The point of contact is the intersection of a horizontal line at \(d\) (baseline) and a linear function with slope \(m\) for the indentation part.

The point of contact is defined as \(\delta=0\) (It’s another fitting parameter).

nanite.poc.poc_fit_constant_polynomial(force, ret_details=False)[source]

Piecewise fit with constant and polynomial

Fit a piecewise function (constant + polynomial) to the baseline and indentation part.

\[F = \frac{\delta^3}{a\delta^2 + b\delta + c} + d\]

This function is defined for all \(\delta>0\). For all \(\delta<0\) the model evaluates to \(d\) (baseline).

I’m not sure where this has been described initially, but it is used e.g. in [RZSK19].

For small indentations, this function exhibits a cubic behavior:

\[y \approx \delta^3/c + d\]

And for large indentations, this function is linear:

\[y \approx \delta/a\]

The point of contact is defined as \(\delta=0\) (It’s another fitting parameter).

nanite.poc.poc_fit_line_polynomial(force, ret_details=False)[source]

Piecewise fit with line and polynomial

Fit a piecewise function (line + polynomial) to the baseline and indentation part.

The linear basline (\(\delta<0\)) is modeled with:

\[F = m \delta + d\]

The indentation part (\(\delta>0\)) is modeled with:

\[F = \frac{\delta^3}{a\delta^2 + b\delta + c} + m \delta + d\]

For small indentations, this function exhibits a linear and only slightly cubic behavior:

\[y \approx \delta^3/c + m \delta + d\]

And for large indentations, this function is linear:

\[y \approx \left( \frac{1}{a} + m \right) \delta\]

The point of contact is defined as \(\delta=0\) (It’s another fitting parameter).

See also

poc_fit_constant_polynomial

polynomial-only version

nanite.poc.poc_frechet_direct_path(force, ret_details=False)[source]

Fréchet distance to direct path

The indentation part is transformed to normalized coordinates (force and corresponding x in range [0, 1]). The point with the largest distance to the line from (0, 0) to (1, 1) is the contact point.

This method is robust with regard to tilted baselines and is a good initial guess for fitting-based POC estimation approaches.

Note that the length of the baseline influences the returned contact point. For shorter baselines, the contact point will be closer to the point of maximum indentation.

nanite.poc.poc_gradient_zero_crossing(force, ret_details=False)[source]

Gradient zero-crossing of indentation part

  1. Apply a moving average filter to the curve

  2. Compute the gradient

  3. Cut off gradient at maximum with a 10 point reserve

  4. Apply a moving average filter to the gradient

  5. The POC is the index of the averaged gradient curve where the values are below 1% of the gradient maximum, measured from the indentation maximum (not from baseline).

nanite.poc.POC_METHODS = [<function poc_deviation_from_baseline>, <function poc_fit_constant_line>, <function poc_fit_constant_polynomial>, <function poc_fit_line_polynomial>, <function poc_frechet_direct_path>, <function poc_gradient_zero_crossing>]

List of all methods available for contact point estimation

Modeling

Methods and constants

nanite.model.compute_anc_parms(idnt, model_key)[source]

Compute ancillary parameters for a force-distance dataset

Ancillary parameters include parameters that:

  • are unrelated to fitting: They may just be important parameters to the user.

  • require the entire dataset: They cannot be extracted during fitting, because they require more than just the approach xor retract curve to compute (e.g. hysteresis, jump of retract curve at maximum indentation). They may, additionally, depend on initial fit parameters set by the user.

  • require a fit: They are dependent on fitting parameters but are not required during fitting.

Notes

If an ancillary parameter name matches that of a fitting parameter, then it is assumed that it can be used for fitting. Please see nanite.indent.Indentation.get_initial_fit_parameters() and nanite.fit.guess_initial_parameters().

Ancillary parameters are set to np.nan if they cannot be computed.

Parameters:
  • idnt (nanite.indent.Indentation) – The force-distance data for which to compute the ancillary parameters

  • model_key (str) – Name of the model

Returns:

ancillaries – key-value dictionary of ancillary parameters

Return type:

collections.OrderedDict

nanite.model.get_anc_parm_keys(model_key)[source]

Return the key names of a model’s ancillary parameters

nanite.model.get_anc_parms(idnt, model_key)[source]
nanite.model.get_init_parms(model_key)[source]

Get initial fit parameters for a model

nanite.model.get_model_by_name(name)[source]

Convenience function to obtain a model by name instead of by key

nanite.model.get_parm_name(model_key, parm_key)[source]

Return parameter label

Parameters:
  • model_key (str) – The model key (e.g. “hertz_cone”)

  • parm_key (str) – The parameter key (e.g. “E”)

Returns:

parm_name – The parameter label (e.g. “Young’s Modulus”)

Return type:

str

nanite.model.get_parm_unit(model_key, parm_key)[source]

Return parameter unit

Parameters:
  • model_key (str) – The model key (e.g. “hertz_cone”)

  • parm_key (str) – The parameter key (e.g. “E”)

Returns:

parm_unit – The parameter unit (e.g. “Pa”)

Return type:

str

Modeling core class

exception nanite.model.core.ModelError[source]
exception nanite.model.core.ModelImplementationError[source]
exception nanite.model.core.ModelImplementationWarning[source]
exception nanite.model.core.ModelImportError[source]
exception nanite.model.core.ModelIncompleteError[source]
class nanite.model.core.NaniteFitModel(model_module)[source]

Initialize the model with an imported Python module

compute_ancillaries(fd)[source]

Compute ancillary parameters for a force-distance dataset

Ancillary parameters include parameters that:

  • are unrelated to fitting: They may just be important parameters to the user.

  • require the entire dataset: They cannot be extracted during fitting, because they require more than just the approach xor retract curve to compute (e.g. hysteresis, jump of retract curve at maximum indentation). They may, additionally, depend on initial fit parameters set by the user.

  • require a fit: They are dependent on fitting parameters but are not required during fitting.

Notes

If an ancillary parameter name matches that of a fitting parameter, then it is assumed that it can be used for fitting. Please see nanite.indent.Indentation.get_initial_fit_parameters() and nanite.fit.guess_initial_parameters().

Ancillary parameters are set to np.nan if they cannot be computed.

Parameters:

fd (nanite.indent.Indentation) – The force-distance data for which to compute the ancillary parameters

Returns:

ancillaries – key-value dictionary of ancillary parameters

Return type:

collections.OrderedDict

get_anc_parm_keys()[source]

Return the key names of a model’s ancillary parameters

get_parm_name(key)[source]

Return parameter label

Parameters:

key (str) – The parameter key (e.g. “E”)

Returns:

parm_name – The parameter label (e.g. “Young’s Modulus”)

Return type:

str

get_parm_unit(key)[source]

Return parameter unit

Parameters:

key (str) – The parameter key (e.g. “E”)

Returns:

parm_unit – The parameter unit (e.g. “Pa”)

Return type:

str

nanite.model.core.compute_anc_max_indent(fd)[source]

Compute ancillary parameter ‘Maximum indentation’

nanite.model.core.ANCILLARY_COMMON = {'max_indent': ('Maximum indentation', 'm', <function compute_anc_max_indent>)}

Common ancillary parameters

Residuals and weighting

nanite.model.residuals.compute_contact_point_weights(cp, delta, weight_dist=5e-07)[source]

Compute contact point weights

Parameters:
  • cp (float) – Fitted contact point value

  • delta (1d ndarray of length N) – The indentation array along which weights will be computed.

  • weight_width (float) – The distance from cp until which weights will be applied.

Returns:

weights – The weights.

Return type:

1d ndarray of length N

Notes

All variables should be given in the same units. The weights increase linearly from increasing distances of delta-cp from 0 to 1 and are 1 outside of the weight width abs(delta-cp)>weight_width.

nanite.model.residuals.get_default_modeling_wrapper(model_function)[source]

Return a wrapper around the default nanite modeling function

nanite.model.residuals.get_default_residuals_wrapper(model_function)[source]

Return a wrapper around the default nanite residual function

nanite.model.residuals.model_direction_agnostic(model_function, params, delta)[source]

Call model_function while making sure that data are in correct order

TODO: Re-evaluate usefulness of this method.

nanite.model.residuals.residual(params, delta, force, model, weight_cp=5e-07)[source]

Compute residuals for fitting

Parameters:
  • params (lmfit.Parameters) – The fitting parameters for model

  • delta (1D ndarray of lenght M) – The indentation distances

  • force (1D ndarray of length M) – The corresponding force data

  • model (callable) – A model function accepting the arguments params and delta

  • weight_cp (positive float or zero/False) – The distance from the contact point until which linear weights will be applied. Set to zero to disable weighting.

nanite.model.weight.weight_cp(*args, **kwargs)[source]

Models

Each model is implemented as a submodule in nanite.model. For instance nanite.model.model_hertz_parabolic. Each of these modules implements the following functions (which are not listed for each model in the subsections below), here with the (non-existent) example module model_submodule:

nanite.model.model_submodule.get_parameter_defaults()

Return the default parameters of the model.

nanite.model.model_submodule.model()

Wrap the actual model for fitting.

nanite.model.model_submodule.residual()

Compute the residuals during fitting (optional).

In addition, each submodule contains the following attributes:

nanite.model.model_submodule.model_doc

The doc-string of the model function.

nanite.model.model_submodule.model_key

The model key used in the command line interface and during scripting.

nanite.model.model_submodule.model_name

The name of the model.

nanite.model.model_submodule.parameter_keys

Parameter keys of the model for higher-level applications.

nanite.model.model_submodule.parameter_names

Parameter names of the model for higher-level applications.

nanite.model.model_submodule.parameter_units

Parameter units for higher-level applications.

Ancillary parameters may also be defined like so:

nanite.model.model_submodule.compute_ancillaries()

Function that returns a dictionary with ancillary parameters computed from an Indentation instance.

nanite.model.model_submodule.parameter_anc_keys

Ancillary parameter keys

nanite.model.model_submodule.parameter_anc_names

Ancillary parameter names

nanite.model.model_submodule.parameter_anc_units

Ancillary parameter units

conical indenter (Hertz)

model key

hertz_cone

model name

conical indenter (Hertz)

model location

nanite.model.model_conical_indenter

nanite.model.model_conical_indenter.hertz_conical(delta, E, alpha, nu, contact_point=0, baseline=0)[source]

Hertz model for a conical indenter

\[F = \frac{2\tan\alpha}{\pi} \frac{E}{1-\nu^2} \delta^2\]
Parameters:
  • delta (1d ndarray) – Indentation [m]

  • E (float) – Young’s modulus [N/m²]

  • alpha (float) – Half cone angle [degrees]

  • nu (float) – Poisson’s ratio

  • contact_point (float) – Indentation offset [m]

  • baseline (float) – Force offset [N]

Returns:

F – Force [N]

Return type:

float

Notes

These approximations are made by the Hertz model:

  • The sample is isotropic.

  • The sample is a linear elastic solid.

  • The sample is extended infinitely in one half space.

  • The indenter is not deformable.

  • There are no additional interactions between sample and indenter.

Additional assumptions:

  • infinitely sharp probe

References

Love (1939) [Lov39], Sneddon (1965) [Sne65] (equation 6.4)

parabolic indenter (Hertz)

model key

hertz_para

model name

parabolic indenter (Hertz)

model location

nanite.model.model_hertz_paraboloidal

nanite.model.model_hertz_paraboloidal.hertz_paraboloidal(delta, E, R, nu, contact_point=0, baseline=0)[source]

Hertz model for a paraboloidal indenter

\[F = \frac{4}{3} \frac{E}{1-\nu^2} \sqrt{R} \delta^{3/2}\]
Parameters:
  • delta (1d ndarray) – Indentation [m]

  • E (float) – Young’s modulus [N/m²]

  • R (float) – Tip radius [m]

  • nu (float) – Poisson’s ratio

  • contact_point (float) – Indentation offset [m]

  • baseline (float) – Force offset [N]

Returns:

F – Force [N]

Return type:

float

Notes

The derivation in [Sne65] reads

\[F = \frac{4}{3} \frac{E}{1-\nu^2} \sqrt{2k} \delta^{3/2},\]

where \(k\) is defined by the paraboloid equation

\[\rho^2 = 4kz.\]

As follows from the derivations in [LL59], this model is valid for either of the two cases:

  • Indentation of a plane (infinite radius) with Young’s modulus \(E\) by a sphere with infinite Young’s modulus and radius \(R\), or

  • Indentation of a sphere with radius \(R\) and Young’s modulus \(E\) by a plane (infinite radius) with infinite Young’s modulus.

These approximations are made by the Hertz model:

  • The sample is isotropic.

  • The sample is a linear elastic solid.

  • The sample is extended infinitely in one half space.

  • The indenter is not deformable.

  • There are no additional interactions between sample and indenter.

Additional assumptions:

  • no surface forces

  • If the indenter is spherical, then its radius \(R\) is much larger than the indentation depth \(\delta\).

References

Sneddon (1965) [Sne65] (equation 6.9), Theory of Elasticity by Landau and Lifshitz (1959) [LL59] (§9 Solid bodies in contact, equation 9.14)

pyramidal indenter, three-sided (Hertz)

model key

hertz_pyr3s

model name

pyramidal indenter, three-sided (Hertz)

model location

nanite.model.model_hertz_three_sided_pyramid

nanite.model.model_hertz_three_sided_pyramid.hertz_three_sided_pyramid(delta, E, alpha, nu, contact_point=0, baseline=0)[source]

Hertz model for three sided pyramidal indenter

\[F = 0.887 \tan\alpha \cdot \frac{E}{1-\nu^2} \delta^2\]
Parameters:
  • delta (1d ndarray) – Indentation [m]

  • E (float) – Young’s modulus [N/m²]

  • alpha (float) – Inclination angle of the pyramidal face [degrees]

  • nu (float) – Poisson’s ratio

  • contact_point (float) – Indentation offset [m]

  • baseline (float) – Force offset [N]

Returns:

F – Force [N]

Return type:

float

Notes

These approximations are made by the Hertz model:

  • The sample is isotropic.

  • The sample is a linear elastic solid.

  • The sample is extended infinitely in one half space.

  • The indenter is not deformable.

  • There are no additional interactions between sample and indenter.

  • The inclination angle of the pyramidal face (in radians) must be close to zero.

References

Bilodeau et al. 1992 [Bil92]

spherical indenter (Sneddon)

model key

sneddon_spher

model name

spherical indenter (Sneddon)

model location

nanite_model_sneddon_spher.model_sneddon_spherical

nanite_model_sneddon_spher.model_sneddon_spherical.delta_of_a(a, R)

Compute indentation from contact area radius (wrapper)

nanite_model_sneddon_spher.model_sneddon_spherical.get_a(R, delta, accuracy=1e-12)

Compute the contact area radius (wrapper)

nanite_model_sneddon_spher.model_sneddon_spherical.hertz_spherical(delta, E, R, nu, contact_point=0.0, baseline=0.0)

Hertz model for Spherical indenter - modified by Sneddon

This model is only available after installing the nanite_model_sneddon_spher Python package.

\[\begin{split}F &= \frac{E}{1-\nu^2} \left( \frac{R^2+a^2}{2} \ln \! \left( \frac{R+a}{R-a}\right) -aR \right)\\ \delta &= \frac{a}{2} \ln \! \left(\frac{R+a}{R-a}\right)\end{split}\]

(\(a\) is the radius of the circular contact area between bead and sample.)

Parameters:
  • delta (1d ndarray) – Indentation [m]

  • E (float) – Young’s modulus [N/m²]

  • R (float) – Tip radius [m]

  • nu (float) – Poisson’s ratio

  • contact_point (float) – Indentation offset [m]

  • baseline (float) – Force offset [N]

Returns:

F – Force [N]

Return type:

float

Notes

These approximations are made by the Hertz model:

  • The sample is isotropic.

  • The sample is a linear elastic solid.

  • The sample is extended infinitely in one half space.

  • The indenter is not deformable.

  • There are no additional interactions between sample and indenter.

Additional assumptions:

  • no surface forces

References

Sneddon (1965) [Sne65] (equations 6.13 and 6.15)

spherical indenter (Sneddon, truncated power series)

model key

sneddon_spher_approx

model name

spherical indenter (Sneddon, truncated power series)

model location

nanite.model.model_sneddon_spherical_approximation

nanite.model.model_sneddon_spherical_approximation.hertz_sneddon_spherical_approx(delta, E, R, nu, contact_point=0, baseline=0)[source]

Hertz model for Spherical indenter - approximation

\[F = \frac{4}{3} \frac{E}{1-\nu^2} \sqrt{R} \delta^{3/2} \left(1 - \frac{1}{10} \frac{\delta}{R} - \frac{1}{840} \left(\frac{\delta}{R}\right)^2 + \frac{11}{15120} \left(\frac{\delta}{R}\right)^3 + \frac{1357}{6652800} \left(\frac{\delta}{R}\right)^4 \right)\]
Parameters:
  • delta (1d ndarray) – Indentation [m]

  • E (float) – Young’s modulus [N/m²]

  • R (float) – Tip radius [m]

  • nu (float) – Poisson’s ratio

  • contact_point (float) – Indentation offset [m]

  • baseline (float) – Force offset [N]

Returns:

F – Force [N]

Return type:

float

Notes

These approximations are made by the Hertz model:

  • The sample is isotropic.

  • The sample is a linear elastic solid.

  • The sample is extended infinitely in one half space.

  • The indenter is not deformable.

  • There are no additional interactions between sample and indenter.

Additional assumptions:

  • no surface forces

Truncated power series approximation:

This model is a truncated power series approximation of spherical indenter (Sneddon). The expected error is more than four magnitues lower than the signal (see e.g. Approximating the Hertzian model with a spherical indenter). The Bio-AFM analysis software by JPK/Bruker uses the same model.

References

Sneddon (1965) [Sne65] (equations 6.13 and 6.15), Dobler (personal communication, 2018) [Dob18]

Fitting

exception nanite.fit.FitDataError[source]
exception nanite.fit.FitKeyError[source]
exception nanite.fit.FitWarning[source]
class nanite.fit.FitProperties[source]

Fit property manager class

Provide convenient access to fit properties as a dictionary and dynamically manage resets due to new initial parameters.

Dynamic properties include:

  • set “params_initial” to None if the “model_key” changes

  • remove all keys except those in FP_DEFAULT if a key that is in FP_DEFAULT changes (All other keys are considered to be obsolete fitting results).

Additional attributes:

reset()[source]
restore(props)[source]

update the dictionary without removing any keys

class nanite.fit.IndentationFitter(idnt, **kwargs)[source]

Fit force-distance curves

Parameters:
  • idnt (nanite.indent.Indentation) – The dataset to fit

  • model_key (str) – A key referring to a model in nanite.model.models_available

  • params_initial (instance of lmfit.Parameters) – Parameters for fitting. If not given, default parameters are used.

  • range_x (tuple of 2) – The range for fitting, see range_type below.

  • range_type (str) –

    One of:

    • absolute:

      Set the absolute fitting range in values given by the x_axis.

    • relative cp:

      In some cases it is desired to be able to fit a model only up until a certain indentation depth (tip position) measured from the contact point. Since the contact point is a fit parameter as well, this requires a two-pass fitting.

  • preprocessing (list of str) – Preprocessing step identifiers

  • preprocessing_options (dict of dicts) – Preprocessing keyword arguments of steps (if applicable)

  • segment (int) – Segment index (e.g. 0 for approach)

  • weight_cp (float) – Weight the contact point region which shows artifacts that are difficult to model with e.g. Hertz.

  • gcf_k (float) – Geometrical correction factor \(k\) for non-single-contact data. The measured indentation is multiplied by this factor to correct for experimental geometries during fitting, e.g. gcf_k=0.5 for parallel-place compression.

  • optimal_fit_edelta (bool) – Search for the optimal fit by varying the maximal indentation depth and determining a plateau in the resulting Young’s modulus (fitting parameter “E”).

  • optimal_fit_num_samples (int) – Number of samples to use for searching the optimal fit

compute_emodulus_vs_mindelta(callback=None)[source]

Compute elastic modulus vs. minimal indentation curve

static compute_opt_mindelta(emoduli, indentations)[source]

Determine the plateau of an emodulus-indentation curve

The following procedure is performed:

  1. Smooth the emodulus data with a Butterworth filter

  2. Label sequences that have similar values by binning into ten regions between the min and max.

  3. Ignore sequences with emodulus that is smaller than the binning size.

  4. Determine the longest sequence.

fit()[source]

Fit the approach-retract data to a model function

get_initial_parameters(idnt=None, model_key='hertz_para')[source]

Get initial fit parameters for a specific model

nanite.fit.guess_initial_parameters(idnt=None, model_key='hertz_para', common_ancillaries=True, model_ancillaries=True)[source]

Guess initial fitting parameters

Parameters:
  • idnt (nanite.indent.Indentation) – The dataset to use for guessing initial fitting parameters using ancillary parameters

  • model_key (str) – The model key

  • common_ancillaries (bool) – Guess global ancillary parameters (such as contact point)

  • model_ancillaries (bool) – Use model-related ancillary parameters

nanite.fit.obj2bytes(obj)[source]

Bytes representation of an object for hashing

Rating

Features

class nanite.rate.features.IndentationFeatures(dataset=None)[source]
static compute_features(idnt, which_type='all', names=None, ret_names=False)[source]

Compute the features for a data set

Parameters:
  • idnt (nanite.Indentation) – A dataset to rate

  • names (list of str) – The names of the rating methods to use, e.g. [“rate_apr_bumps”, “rate_apr_mon_incr”]. If None (default), all available rating methods are used.

Notes

names may include features that are excluded by which_type. E.g. if a “bool” feature is in names but which_type is “float”, then the “bool” feature will be silently ignored.

feat_bin_apr_spikes_count()[source]

spikes during IDT

Sudden spikes in indentation curve

feat_bin_cp_position()[source]

CP outside of data range

Contact point position outside of range

feat_bin_size()[source]

dataset too small

Number of points in indentation curve

feat_con_apr_flatness()[source]

flatness of APR residuals

fraction of the positive-gradient residuals in the approach part

feat_con_apr_size()[source]

relative APR size

length of the approach part relative to the indentation part

feat_con_apr_sum()[source]

residuals of APR

absolute sum of the residuals in the approach part

feat_con_bln_slope()[source]

slope of BLN

slope obtained from a linear least-squares fit to the baseline

feat_con_bln_variation()[source]

variation in BLN

comparison of the forces at the beginning and at the end of the baseline

feat_con_cp_curvature()[source]

curvature at CP

curvature of the force-distance data at the contact point

feat_con_cp_magnitude()[source]

residuals at CP

mean value of the residuals around the contact point

feat_con_idt_maxima_75perc()[source]

maxima in IDT residuals

sum of the indentation residuals’ maxima in three intervals in-between 25% and 100% relative to the maximum indentation

feat_con_idt_monotony()[source]

monotony of IDT

change of the gradient in the indentation part

feat_con_idt_spike_area()[source]

area of IDT spikes

area of spikes appearing in the indentation part

feat_con_idt_sum()[source]

overall IDT residuals

sum of the residuals in the indentation part

feat_con_idt_sum_75perc()[source]

residuals in 75% IDT

sum of the residuals in the indentation part in-between 25% and 100% relative to the maximum indentation

classmethod get_feature_funcs(which_type='all', names=None)[source]

Return functions that compute features from a dataset

Parameters:
  • names (list of str) – The names of the rating methods to use, e.g. [“rate_apr_bumps”, “rate_apr_mon_incr”]. If None (default), all available rating methods are returned.

  • which_type (str) – Which features to return: [“all”, “bool”, “float”].

Returns:

raters – Each item in the list consists contains the name of the rating method and the corresponding rating method.

Return type:

list of tuples (name, callable)

classmethod get_feature_names(which_type='all', names=None, ret_indices=False)[source]

Get features names

Parameters:
  • which_type (str or list of str) – Return only features that are of a certain type. See VALID_FEATURE_TYPES for valid strings.

  • names (list of str) – Return only features that are in this list.

  • ret_indices (bool) – If True, also return the internal feature indices.

Returns:

name_list – List of feature names (callables of this class)

Return type:

list of str

property contact_point
property datafit_apr
property datares_apr
dataset

current dataset from which features are computed

property datax_apr
property datay_apr
property has_contact_point
property is_fitted
property is_valid
property meta
nanite.rate.features.VALID_FEATURE_TYPES = ['all', 'binary', 'continuous']

Valid keyword arguments for feature types

Rater

class nanite.rate.rater.IndentationRater(regressor=None, scale=None, lda=None, training_set=None, names=None, weight=True, sample_weight=None, *args, **kwargs)[source]

Rate quality

Parameters:
  • regressor (sciki-learn RegressorMixin) – The regressor used for rating

  • scale (bool) – If True, apply a Standard Scaler. If a regressor based on decision trees is used, the Standard Scaler is not used by default, otherwise it is.

  • lda (bool) – If True, apply a Linear Discriminant Analysis (LDA). If a regressor based on a decision tree is used, LDA is not used by default, otherwise it is.

  • training_set (tuple of (X, y)) – The training set (samples, response)

  • names (list of str) – Feature names to use

  • weight (bool) – Weight the input samples by the number of occurrences or with sample_weight. For tree-based classifiers, set this to True to avoid bias.

  • sample_weight (list-like) – The sample weights. If set to None sample weights are computed from the training set.

  • *args (list) – Positional arguments for IndentationFeatures

  • **kwargs – Keyword arguments for IndentationFeatures

See also

sklearn.preprocessing.StandardScaler

Standard scaler

sklearn.discriminant_analysis.LinearDiscriminantAnalysis

Linear discriminant analysis

nanite.rate.regressors.reg_trees

List of regressors that are identified as tree-based

static compute_sample_weight(X, y)[source]

Weight samples according to occurrence in y

static get_training_set_path(label='zef18')[source]

Return the path to a training set shipped with nanite

Training sets are stored in the nanite.rate module path with ts_ prepended to label.

classmethod load_training_set(path: Path | str = None, names: List[str] = None, which_type: Literal['all', 'binary', 'continuous'] | List = None, replace_inf: bool = True, impute_zero_rated_nan: bool = True, remove_nan: bool = True, ret_names: bool = False)[source]

Load a training set from a directory

Parameters:
  • path (pathlib.Path or str) – Optional path to the training set directory. If none is specified, the default “zef18” is loaded.

  • names (list of str) – List of features to use, defaults to all features.

  • which_type (str) – Which type of feature to return see VALID_FEATURE_TYPES for valid options. By default, only the “continuous” features are imported. The “binary” features are not needed for training; they are used to sort out new force-distance data.

  • replace_inf (bool) – Replace infinity-valued feature values with 2 * sign * max(abs(values)).

  • impute_zero_rated_nan (bool) – If there are nan-valued features that have a zero response (rated worst), replace those feature values with the mean of the zero-response features that are not nan-valued.

  • remove_nan (bool) – Remove any nan-valued features (after impute_zero_rated_nan was applied). This is necessary, since skimage cannot handle nan-valued sample values.

  • ret_names (bool) – Return the names of the features in addition to the samples and response.

Returns:

  • samples (2d ndarray) – Sample values with axes (data_size, num_features)

  • response (1d ndarray) – Response array of length data_size

  • names (list, optional) – List of feature names corresponsing to axis 1 in samples

rate(samples=None, datasets=None)[source]

Perform rating step

Parameters:
  • samples (1d or 2d ndarray (cast to 2d ndarray) or None) – Measured samples, if set to None, dataset must be given.

  • dataset (list of nanite.Indentation) – Full, fitted measurement

Returns:

ratings – Resulting ratings

Return type:

list

names

feature names used by the regressor pipeline

pipeline

sklearn pipeline with transforms (and regressor if given)

nanite.rate.rater.get_available_training_sets()[source]

List of internal training sets

nanite.rate.rater.get_rater(regressor, training_set='zef18', names=None, lda=None, **reg_kwargs)[source]

Convenience method to get a rater

Parameters:
  • regressor (str or RegressorMixin) – If a string, must be in reg_names.

  • training_set (str or pathlib.Path or tuple (X, y)) – A string label representing a training set shipped with nanite, the path to a training set, or a tuple representing the training set (samples, response) for use with sklearn.

  • names (list of str) – Only use these features for rating

  • lda (bool) – Perform linear discriminant analysis

Returns:

irater – The rating instance.

Return type:

nanite.IndentationRater

Regressors

scikit-learn regressors and their keyword arguments

nanite.rate.regressors.reg_names = ['AdaBoost', 'Decision Tree', 'Extra Trees', 'Gradient Tree Boosting', 'Random Forest', 'SVR (RBF kernel)', 'SVR (linear kernel)']

List of available default regressor names

nanite.rate.regressors.reg_trees = ['AdaBoostRegressor', 'DecisionTreeRegressor', 'ExtraTreesRegressor', 'GradientBoostingRegressor', 'RandomForestRegressor']

List of tree-based regressor class names (used for keyword defaults in IndentationRater)

Manager

Save and load user-rated datasets

class nanite.rate.io.RateManager(path, verbose=0)[source]

Manage user-defined rates

export_training_set(path)[source]
get_cross_validation_score(regressor, training_set=None, n_splits=20, random_state=42)[source]

Regressor cross-validation scoring

Cross-validation is used to identify regressors that over-fit the train set by splitting the train set into multiple learn/test sets and quantifying the regressor performance for each split.

Parameters:
  • regressor (str or RegressorMixin) – If a string, must be in reg_names.

  • training_set (X, y) – If given, do not use self.samples

Notes

A skimage.model_selection.KFold cross validator is used in combination with the mean squared error score.

Cross-validation score is computed from samples that are filtered with the binary features and only from samples that do not contain any nan values.

get_rates(which='user', training_set='zef18')[source]
which: str

Which rating to return: “user” or a regressor name

get_training_set(which_type='all', prefilter_binary=False, remove_nans=False, transform=False)[source]

Return (X, y) training set

property datasets
path

Path to the manual ratings (directory or .h5 file)

property ratings
property samples

The individual sample ratings computed by afmlib

verbose

verbosity level

nanite.rate.io.hash_file(path, blocksize=65536)[source]

Compute sha256 hex-hash of a file

Parameters:
  • path (str or pathlib.Path) – path to the file

  • blocksize (int) – block size read from the file

Returns:

hex – The first six characters of the hash

Return type:

str

nanite.rate.io.hdf5_rated(h5path, indent)[source]

Test whether an indentation has already been rated

Return type:

is_rated, rating, comment

nanite.rate.io.load(path, meta_only=False, verbose=0)[source]

Notes

The .fit_properties attribute of each Indentation instance is overridden by a simple dictionary, so its functionalities are not available anymore.

nanite.rate.io.load_hdf5(path, meta_only=False)[source]
nanite.rate.io.save_hdf5(h5path, indent, user_rate, user_name, user_comment, h5mode='a')[source]

Store all relevant data of a user rating into an hdf5 file

Parameters:
  • h5path (str or pathlib.Path) – Path to HDF5 rating container where data will be stored

  • indent (nanite.Indentation) – The experimental data processed and fitted with nanite

  • user_rate (float) – Rating given by the user

  • user_name (str) – Name of the rating user

Quantitative maps

exception nanite.qmap.DataMissingWarning[source]
class nanite.qmap.QMap(path_or_group, meta_override=None, callback=None)[source]

Quantitative force spectroscopy map handling

Parameters:
  • path_or_group (str or pathlib.Path or afmformats.afm_group.AFMGroup) – The path to the data file or an instance of AFMGroup

  • meta_override (dict) – Dictionary with metadata that is used when loading the data in path.

  • callback (callable or None) – A method that accepts a float between 0 and 1 to externally track the process of loading the data.

static feat_fit_contact_point(idnt)[source]

Contact point of the fit

static feat_fit_youngs_modulus(idnt)[source]

Young’s modulus

static feat_meta_rating(idnt)[source]

Rating