Implement and document the new Wöhler analysis strategy

CHANGELOG.md

@@ -9,6 +9,7 @@ In this file noteworthy changes of new releases of pyLife are documented since
 ### Breaking changes
 
 * Drop support for python-3.8 ast it is end of life
+* The Wöhler analysis slightly changed (#138)
 
 
 ### Minor improvements

docs/materialdata/woehler.rst

@@ -1,6 +1,19 @@
 The ``woehler`` module
 ######################
 
+Notes on Wöhler data analysis
+=============================
+
+Even though there are established Wöhler analysis algorithms they all leave
+some degrees of freedom.  Therfore we document the way we perform the analysis
+in pyLife in the following section.
+
+.. toctree::
+   :maxdepth: 1
+
+   woehler/notes_on_analysis
+
+
 Module description
 ==================
 
@@ -29,7 +42,6 @@ Analyzers
    woehler/probit
    woehler/maxlikeinf
    woehler/maxlikefull
-   woehler/bayesian
 
 
 Helpers

docs/materialdata/woehler/notes_on_analysis.rst

@@ -0,0 +1,110 @@
+Notes on Wöhler analysis
+########################
+
+Even though there are established Wöhler analysis algorithms they all leave
+some degrees of freedom.  Therfore we document the way we perform the analysis
+in pyLife in this section.
+
+
+Basic principles of analysis
+============================
+
+There are some minimal requirements for the data to perform a Wöhler analysis
+in a meaningful way. Some of them are mandatory and we outright refuse the
+analysis if there are not met. Some we do not mandate but warn, because it is
+likely that the result is not meaningful and should not be used for any serious
+purpose.
+
+
+The concept of load levels
+--------------------------
+
+In general, fatigue data is the cyclic load plotted over the cycles endured.
+In order to perform an analysis, we need what we call load levels.  So the load
+values chosen for the fatigue data is not arbitrary, there should always be
+multiple tests – ideally five – on a single exact same load value. Those we
+call the load levels.
+
+
+The finite and infinite zone
+----------------------------
+
+We separate the Wöhler data in two zones. The finite zone is all the load
+levels at which all specimens fracture. The infinite zone is all the load
+levels at which there is at least one runout. Usually in order to perform a
+proper Wöhler analysis you need at least two load levels in each zone. For the
+infinite zone that means at least two mixed levels that has at least one runout
+and at least one fracture.
+
+The analysis of Wöhler data basically works in two steps. First we calculate
+the finite part, that is the sloped part with the fractures to calculate the
+slope `k` the scatter in life time direction `TN` and the imaginary 50% failure
+load at one cycle. Then we calculate the infinite values, those are the
+endurance limit in load direction `SD` and the scatter in load direction `SD`.
+
+
+Calculating the finite values
+-----------------------------
+
+Calculating the slope `k` is basically a linear regression of all the fractures
+in the finite zone.
+
+
+Steps of the analysis
+=====================
+
+The analysis of Wöhler data basically works in two steps. First we calculate
+the finite part, that is the sloped part with the fractures to calculate the
+slope `k` the scatter in life time direction `TN` and the imaginary 50% failure
+load at one cycle. Then we calculate the infinite values, those are the
+endurance limit in load direction `SD` and the scatter in load direction `SD`.
+
+
+Calculating the finite values
+-----------------------------
+
+Calculating the slope `k` is in the easiest case a linear regression of all the
+fractures in the finite zone. That's what is done by
+:class:`~pylife.materialdata.woehler.Elementary`.
+
+
+Calculating the infinite values
+-------------------------------
+
+The infinite values cannot be calculated by a linear regression. There are two
+established methods to calculate them:
+
+* the maximum likelihood method implemented in
+  :class:`~pylife.materialdata.woehler.MaxLikeFull` and
+  :class:`~pylife.materialdata.woehler.MaxLikeInf`.
+* the probit method implemented in :class:`~pylife.materialdata.woehler.Probit`
+
+Please refer to the documentation of the implementing modules for details.
+
+
+Likelihood estimation strategy
+==============================
+
+The maximum likelihood methods basically try to find a solution that maximizes
+the likelihood that the solution is correct.  So it needs to a way of
+estimating the likelihood for a solution in question.  Again we do that
+separately for the finite and the infinite values.  That raises the question,
+which load levels to take into account to evaluate the likelihood of the finite
+resp. infinite zone.  We have three options implemented for that, all of which
+use all fractures to evaluate the infinite zone.  For the finite zone they do
+the following
+
+1. By default we use all fractures on pure fracture levels and the fractures of
+   the highest mixed level.
+   :class:`~pylife.materialdata.woehler.likelihood.LikelihoodHighestMixedLevel`
+
+2. The second option is to use only the fractures of the pure fracture levels.
+   :class:`~pylife.materialdata.woehler.likelihood.LikelihoodPureFiniteZone`
+
+3. The third option is to use all fractures.
+   :class:`~pylife.materialdata.woehler.likelihood.LikelihoodAllFractures`
+
+
+You can also implement your own likelihood estimator by subclassing
+:class:`~pylife.materialdata.woehler.likelihood.AbstractLikelihood` and
+then

src/pylife/materialdata/woehler/elementary.py

@@ -18,12 +18,18 @@
 import pandas as pd
 import scipy.stats as stats
 
-from .likelihood import Likelihood
+from .likelihood import (
+    LikelihoodAllFractures,
+    LikelihoodHighestMixedLevel,
+    LikelihoodPureFiniteZone,
+    LikelihoodLegacy
+)
 from .pearl_chain import PearlChainProbability
 import pylife.utils.functions as functions
 from . import FatigueData, determine_fractures
 import warnings
 
+
 class Elementary:
     """Base class to analyze SN-data.
 
@@ -46,10 +52,68 @@ def __init__(self, fatigue_data):
            The SN-data to be analyzed.
         """
         self._fd = self._get_fatigue_data(fatigue_data)
-        self._lh = self._get_likelihood()
+        self.use_highest_mixed_level()
+
+    def use_old_likelihood_estimation(self):
+        """Use the old (until pyLife-2.1.x) likelihood estimation.
+
+        That uses all fractures for the finite likelihood and only the mixed levels for
+        the infinite.
+
+        Returns
+        -------
+        self
+        """
+        self._lh = LikelihoodLegacy(self._fd)
+        return self
+
+    def use_highest_mixed_level(self):
+        """Use fractures of pure fracture levels and the highest mixed level
+        to determine the likelihood in the finite zone.
 
-    def _get_likelihood(self):
-        return Likelihood(self._fd)
+        This is the default
+
+        Returns
+        -------
+        self
+        """
+        self._lh = LikelihoodHighestMixedLevel(self._fd)
+        return self
+
+    def use_all_fractures(self):
+        """Use all fractures to determine the likelihood in the finite zone.
+
+        Returns
+        -------
+        self
+        """
+        self._lh = LikelihoodAllFractures(self._fd)
+        return self
+
+    def use_only_pure_fracture_levels(self):
+        """Use only fractures of pure fracture levels to determine the likelihood in the finite zone.
+
+        Returns
+        -------
+        self
+        """
+        self._lh = LikelihoodPureFiniteZone(self._fd)
+        return self
+
+    def use_custom_likelihood_estimation(self, likelihood_class):
+        """Inject a custom Likelihood calculation class to determine the likelihood.
+
+        Parameters
+        ----------
+        likelihood : Class implementing :class:`~pylife.materialdata.woehler.likelihood.Likelihood`
+           The likelihood calculation class
+
+        Returns
+        -------
+        self
+        """
+        self._lh = likelihood_class(self._fd)
+        return self
 
     def _get_fatigue_data(self, fatigue_data):
         if isinstance(fatigue_data, pd.DataFrame):
@@ -138,11 +202,14 @@ def _specific_analysis(self, wc):
     def bayesian_information_criterion(self):
         """The Bayesian Information Criterion
 
-        Bayesian Information Criterion is a criterion for model selection among
+        Bayesian Information Criterion (BIC) is a criterion for model selection among
         a finite set of models; the model with the lowest BIC is preferred.
         https://www.statisticshowto.datasciencecentral.com/bayesian-information-criterion/
 
         Basically the lower the better the fit.
+
+        Note that the BIC is not suitable to compare the results of two different
+        likelihoods.
         """
         if not hasattr(self,"_bic"):
             raise ValueError("BIC value undefined. Analysis has not been conducted.")

src/pylife/materialdata/woehler/fatigue_data.py

@@ -80,6 +80,10 @@ def cycles(self):
         '''the cycle numbers'''
         return self._obj.cycles
 
+    @property
+    def fracture(self):
+        return self._obj.fracture
+
     @property
     def finite_infinite_transition(self):
         '''The start value of the load endurance limit.
@@ -209,6 +213,7 @@ def _calc_finite_zone_manual(self, limit):
         self._finite_zone = self.fractures[self.fractures.load > limit]
         self._infinite_zone = self._obj[self._obj.load <= limit]
 
+
 def determine_fractures(df, load_cycle_limit=None):
     '''Adds a fracture column according to defined load cycle limit