13288720301

Committed 12 Feb 2025 03:19PM UTC coverage: 67.793% (-9.5%) from 77.332%

Build # 13288720301

Build Type

Pull #330

github

Committed by

web-flow

Commit Message

Merge ca74a0d3b into 7e5e3f4f2

Pull Request Pull Request #330: Pr 316

Run Details

38 of 196 new or added lines in 7 files covered. (19.39%)

3 existing lines in 3 files now uncovered.

3753 of 5536 relevant lines covered (67.79%)

6.1 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

16.67

/src/spotpy/algorithms/efast.py

# -*- coding: utf-8 -*-
"""
EFAST algorithm transfered from CRAN R package 'fast' by Dominik Reusser (2015)
Copy of most functions and descriptions and implementation in SPOTPY framework by Anna Herzog (2024)

This Version of FAST differs from the previously impemented Version, as it does not require a precalculation of minimum model runs by the user. It rather uses a hard coded definition of minimum runs based on Cukier or McRae.
It therefore requieres considreably less model runs (for 5 Parameters 2245 (fast) vs. 71 runs (efast)).

For further information on the functions, the origninal R package and the fast algorithm please have a look at the R package description from CRAN, or see the following references:

Reusser, Dominik E., Wouter Buytaert, and Erwin Zehe. 
"Temporal dynamics of model parameter sensitivity for computationally expensive models with FAST (Fourier Amplitude Sensitivity Test)." Water Resources Research 47 (2011): W07551.

Further References:
CUKIER, R. I.; LEVINE, H. B. & SHULER, K. E. Non-Linear 
"Sensitivity Analysis Of Multi-Parameter Model Systems" Journal Of Computational Physics, 1978 , 26 , 1-42
CUKIER, R. I.; FORTUIN, C. M.; SHULER, K. E.; PETSCHEK, A. G. & SCHAIBLY, J. H. 
"Study Of Sensitivity Of Coupled Reaction Systems To Uncertainties In Rate Coefficients .1. Theory" Journal Of Chemical Physics, 1973 , 59 , 3873-3878
SCHAIBLY, J. H. & SHULER, K. E. 
"Study Of Sensitivity Of Coupled Reaction Systems To Uncertainties In Rate Coefficients .2. Applications" Journal Of Chemical Physics, 1973 , 59 , 3879-3888
CUKIER, R. I.; SCHAIBLY, J. H. & SHULER, K. E. 
"Study Of Sensitivity Of Coupled Reaction Systems To Uncertainties In Rate Coefficients .3. Analysis Of Approximations" Journal Of Chemical Physics, 1975 , 63 , 1140-1149


Copyright (c) 2018 by Tobias Houska
This file is part of Statistical Parameter Optimization Tool for Python(SPOTPY).
:author: Anna Herzog, Tobias Houska, Dominik Reusser
"""

import warnings

import numpy as np

from ..analyser import efast_sensitivity, get_modelruns
from . import _algorithm


class efast(_algorithm):
    """
    EFAST Algorithm for (distributed) parameter Sensitivity after FAST algorithm according to Cukier 1975 or McRae 1982
    The Fourier Amplitude Sensitivity Test (FAST) is a method to determine global sensitvities of a model on parameter changes
    within relatively few model runs. 
    """

    _unaccepted_parameter_types = ()

    def __init__(self, *args, **kwargs):
        """
        Input
        ----------
        spot_setup: class
            model: function
                Should be callable with a parameter combination of the parameter-function
                and return an list of simulation results (as long as evaluation list)
            parameter: function
                When called, it should return a random parameter combination. Which can
                be e.g. uniform or Gaussian
            objectivefunction: function
                Should return the objectivefunction for a given list of a model simulation and
                observation.
            evaluation: function
                Should return the true values as return by the model.

        dbname: str
            * Name of the database where parameter, objectivefunction value and simulation results will be saved.

        dbformat: str
            * ram: fast suited for short sampling time. no file will be created and results are saved in an array.
            * csv: A csv file will be created, which you can import afterwards.

        parallel: str
            * seq: Sequentiel sampling (default): Normal iterations on one core of your cpu.
            * mpi: Message Passing Interface: Parallel computing on cluster pcs (recommended for unix os).

        save_sim: boolean
            * True:  Simulation results will be saved
            * False: Simulation results will not be saved
        """

        

        kwargs["algorithm_name"] = "EFast Sampler"
        super(efast, self).__init__(*args, **kwargs)

    # hard coded repetition and frequency values from R package Fast
    d_m_cukier75 = [4,8,6,10,20,22,32,40,38,26,56,62,46,76,96,60,86,126,134,112,92,128,154,196,34,416,106,208,328,198,382,88,348,186,140,170,284,568,302,438,410,248,448,388,596,217,100,488,166]
    min_runs_cukier75 = [np.nan, np.nan, 19, 39, 71, 91, 167, 243, 315, 403, 487, 579, 687, 907, 1019, 1223, 1367, 1655, 1919, 2087, 2351, 2771, 3087, 3427, 3555, 4091, 4467, 4795, 5679, 5763, 6507, 7103, 7523, 8351, 9187, 9667, 10211, 10775, 11339, 7467, 12891, 13739, 14743, 15743, 16975, 18275, 18927, 19907, 20759, 21803]
    omega_m_cukier75 = [np.nan, np.nan, 1, 5, 11, 1, 17, 23, 19, 25, 41, 31, 23, 87, 67, 73, 58, 143, 149, 99, 119, 237, 267, 283, 151, 385, 157, 215, 449, 163, 337, 253, 375, 441, 673, 773, 875, 873, 587, 849, 623, 637, 891, 943, 1171, 1225, 1335, 1725, 1663, 2019]
    d_m_mcrae82 = [4, 8, 6, 10, 20, 22, 32, 40, 38, 26, 56, 62, 46, 76, 96]
    min_runs_mcrae82 = [0, 15, 27, 47, 79, 99, 175, 251, 323, 411, 495, 587, 695, 915, 1027]
    omega_m_mcrae82 = [0, 3, 1, 5, 11, 1, 17, 23, 19, 25, 41, 31, 23, 87, 67]

    def freq_cukier(self, m, i = 1, omega_before = -1):
    
        """
        This function creates an array of independant frequencies accroding to 
        Cukier 1975 for usage in the eFAST method.

        Parameters
        ----------
        m: int
            number of parameters (frequencies) needed
        i: (intern) int
            internally used recursion counter
        omega_before: (intern) int
            internally used previous frequency

        Raises
        ----------
        Exeption:
            "Parameter number m to high, not implemented"
            Execution is stopped, if the number of Parameters is to high (max. Number of Parameters for Cukier = 50)    

        Returns
        ----------
        np.array of float:
            A 1d-Array of independant frequencies to the order of 4

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog

        """

        if i <= 1:
            if m >= len(self.omega_m_cukier75):
                raise Exception("Parameter number m to high, not implemented") 
            else:
                o = self.omega_m_cukier75[m-1]
                return np.append(o, self.freq_cukier(m, i+1, o))
        else:
            o = omega_before + self.d_m_cukier75[m-i]
            if i == m:
                return o
            else:
                return np.append(o, self.freq_cukier(m, i+1, o))
        

    def freq_mcrae82(self, m, i = 1, omega_before = -1):
    
        """
        This function creates an array of independant frequencies accroding to 
        McRae 1982 for usage in the eFAST method.

        Parameters
        ----------
        m: int
            number of parameters (frequencies) needed
        i: (intern) int
            internally used recursion counter
        omega_before: (intern) int
            internally used previous frequency

        Raises
        ----------
        Exeption:
            "Parameter number m to high, not implemented"
            Execution is stopped, if the number of Parameters is to high (max. Number of Parameters for McRae = 15)    

        Returns
        ----------
        np.array of float:
            A 1d-Array of independant frequencies to the order of 4

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog

        """

        if i <= 1:
            if m >= len(self.omega_m_mcrae82):
                raise Exception("Parameter number m to high, not implemented") 
            else:
                o = self.omega_m_mcrae82[m-1]
                return np.append(o, self.freq_mcrae82(m, i+1, o))
        else:
            o = omega_before + self.d_m_mcrae82[m-i]
            if i == m:
                return o
            else:
                return np.append(o, self.freq_mcrae82(m, i+1, o))
        

    def s(self, m, freq = 'cukier'):

        """
        Function that generates a number of equally spaced values between -p1/2 
        and pi/2. The number is determined by the number of runs required of the
        eFAST method for a number of parameters (min_runs_cukier or min_runs_mcrae)

        Parameters
        ----------
        m: int
            number of parameters/ frequencies

        freq: (optional) str
            indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
            Default is Cukier

        Raises
        ----------
        Exeption:
            "Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!"
            Execution is stopped if an invalid method for the frequency selection is provided

        Returns
        ----------
        2D array of float:
            an array of equally spaced values between -pi/2 and pi/2

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog

        """

        if freq == 'cukier':
            min_runs = self.min_runs_cukier75
        elif freq == 'mcrae':
            min_runs = self.min_runs_mcrae82
        else:
            raise Exception("Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!")
        
        r = min_runs[m-1]
        r_range = np.array(range(1,r+1))
        s = np.pi/r * (2*r_range-r-1)/2
    
        return s


    def S(self, m, freq = 'cukier'):
    
        """
        Function to generate an array of values, which provides the base for the parameterdistribution 
        in the FAST method. It is usally not used directly but called from the 
        fast_parameters function.

        Parameters
        ----------
        m: int
            number of parameters/frequencies

        freq: (optional) str
            indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
            Default is Cukier

        Raises
        ----------
        Exeption:
            "Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!"
            Execution is stopped if an invalid method for the frequency selection is provided

        Returns
        ----------
        2D array of float:
            an array with the shape (number of runs, parameters)

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog
        """

        if freq == 'cukier':
            omega = self.freq_cukier(m)
        elif freq == 'mcrae':
            omega = self.freq_mcrae82(m)
        else:
            raise Exception("Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!")
    
        tab = np.outer(self.s(m, freq), omega)
    
        toreturn = np.arcsin(np.sin(tab))/np.pi
    
        # naming array dimensions is not possible with numpy but convention would be toreturn.shape () = (runs, parameternames)
    
        return toreturn


    def rerange(self, data, min_goal = 0, max_goal = 1, center = np.nan):

        """
        This function performes a linear transformation of the data, such that
        afterwards range(data) = (theMin, theMax)

        Parameters
        ----------
        data: array of flaot
            an array with the data to transform
            in this case the parameter distribution generated by function S
        min_goal: float
            the new minimum value (lower parameter bound)
        max_goal: float
            the new maximum value (upper parameter bound)
        center: (optional) float
            indicates which old value should become the new center
            default: (max_goal+min_goal)/2

        Returns
        ----------
        2D array of float:
            array with the transformed data

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog
    
        """
    
        min_data = min(data)
        max_data = max(data)
    
        if np.isnan(center):
            max_data = max_data
            dat = data - min_data
            dat = dat/(max_data-min_data) * (max_goal-min_goal)
            dat = dat + min_goal
            return(dat)
        else:
            # split linear transformation, the data is split into two segments, one blow center, one above, 
            # each segment undergoes it's own linear transformation
            below = data <= center
            above = data > center
            new_data = np.copy(data)
            np.place(new_data, below, self.rerange(np.insert(data[below], 0, center), min_goal = min_goal, max_goal= max_goal/2)[1:])
            np.place(new_data, above, self.rerange(np.insert(data[above], 0, center), min_goal = max_goal/2, max_goal= max_goal)[1:])
            return(new_data)
    

    def fast_parameters(self, minimum, maximum, freq = 'cukier', logscale = np.nan, names = np.nan):
    
        """
        Function for the creation of a FAST Parameter set based on thier range

        Parameters
        ----------
        minimum: array of float 
            array containing the lower parameter boundries
        maximum: array of float 
            array containing the upper parameter boundries
        names: (optional) str
            array containing the parameter names
        freq: (optional) str
            indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
            Default is Cukier
        logscale: (optional) bool
            array containing bool values indicating weather a parameter is varied
            on a logarithmic scale. In that case minimum and maximum are exponents

        Raises
        ----------
        Exeption:
            "Expecting minimum and maximum of same size"
            Execution is stopped if the number of minimum, maximum or logscale values differ

        Returns
        ----------
        value: array of float
            an array with the shape (number of runs, parameters) containing the 
            required parameter sets

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog

        """
        
        if np.isnan(logscale):
            logscale = np.full(minimum.shape, False)
        if np.isnan(names):
            names = ["P"+str(i) for i in range(minimum.shape[0])]
            names = np.array(names)  
    
        n = len(minimum)
    
        if (n != len(maximum)):
            raise Exception("Expecting minimum and maximum of same size") 
        elif(n != len(names)):
            raise Exception("Expecting minimum and names of same size") 
        elif(n != len(logscale)):
            raise Exception("Expecting minimum and logscale of same size") 
    
        toreturn = self.S(m=n, freq = freq) #par_names = names, 
    
        for i in range(0,n):
            toreturn[:,i] = self.rerange(toreturn[:,i], minimum[i], maximum[i])
            if logscale[i]:
                toreturn[:,i] = 10**toreturn[:,i]
            
        return toreturn


    def sample(self, repetitions, freq = 'cukier', logscale = np.nan):
        """
        Samples from the eFAST algorithm.

        Input
        ----------
        repetitions: int
            Maximum number of runs.
        freq: (optional) str
            indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
            Default is Cukier
        logscale: (optional) bool
            array containing bool values indicating weather a parameter is varied
            on a logarithmic scale. In that case minimum and maximum are exponents

        Raises
        ----------
        Warning:
            "specified number of repetitions is smaller than minimum required for eFAST analysis!
            Simultions will be perfomed but eFAST analysis might not be possible"

            If the specified number of repetitions is smaller than required by the algorithm, the
            algortihm will still execute, but the usage of spotpy.analyser.efast.calc_sensitivity() 
            might not be possible. Rather define a large number of reps as only required minimum runs
            will be performed anyway. 
        
        Warning:
            "Number of specified repetitions equals or exeeds number of minimum required repetitions for eFAST analysis
             programm will stop after required runs."

        Returns
        ----------
        database:
            spotpy data base in chosen format with objective function values, parameter values and simulation results

        Author
        ----------
        Tobias Houska, Anna Herzog
        """

        print("generating eFAST Parameters")
        # Get the names of the parameters to analyse
        names = self.parameter()["name"]
        # Get the minimum and maximum value for each parameter from the
        # distribution
        parmin, parmax = self.parameter()["minbound"], self.parameter()["maxbound"]
        self.numberf = len(parmin)
        min_runs = self.min_runs_cukier75[self.numberf-1]

        if min_runs > repetitions:
            warnings.warn("specified number of repetitions is smaller than minimum required for eFAST analysis!\n"+
                          "Simultions will be perfomed but eFAST analysis might not be possible")
        else:
            repetitions = min_runs
            print("Number of specified repetitions equals or exeeds number of minimum required repetitions for eFAST analysis\n"+
                  "programm will stop after required runs.")

        self.set_repetiton(repetitions)

        # Generate Matrix with eFAST parameter sets
        N = self.fast_parameters(parmin, parmax, freq, logscale)

        print("Starting the eFast algorithm with {} repetitions...".format(repetitions))

        # A generator that produces parametersets if called
        param_generator = (
            (rep, N[rep,:]) for rep in range(int(repetitions))
        )
        for rep, randompar, simulations in self.repeat(param_generator):
            # A function that calculates the fitness of the run and the manages the database
            self.postprocessing(rep, randompar, simulations)

            if self.breakpoint == "write" or self.breakpoint == "readandwrite":
                if rep >= lastbackup + self.backup_every_rep:
                    work = (rep, N[rep,:])
                    self.write_breakdata(self.dbname, work)
                    lastbackup = rep
        self.final_call()

    def calc_sensitivity(self, results, dbname, freq = 'cukier'):
        """
        Function to calcultae the sensitivity for a data series (eg. a time series) 
        based on the eFAST algorithm.

        Parameters
        ----------
        results: np.array of void
            spopty restults array as loaded with spotpy.analyser.load_csv_results()
        dbname: (optional) str
            name of file to store sensitivity values        
        freq: (optional) str
            indicates weather to use the frequencies after Cukier 'cukier' or McRae 'mcrae'
            Default is Cukier

        Raises
        ----------
        Exeption:
            "Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!"
            Execution is stopped if an invalid method for the frequency selection is provided

        Returns
        ----------
        database:
            database containing the temporal parameter sensitivities with seperate column
            for each parameter and row for eg. time step

        Author
        ----------
        Dominic Reusser
        spotpy translation: Anna Herzog
    
        """ 
        numberf = len(self.parameter()["minbound"])

        # sens_data = np.full((len(results[0]), numberf), np.nan)
        
        # idendtify frequency
        if freq == 'cukier':
            t_runs = self.min_runs_cukier75[numberf-1]
            t_freq = self.freq_cukier(numberf)
        elif freq == 'mcrae': 
            t_runs = self.min_runs_mcrae82[numberf-1]
            t_freq = self.freq_mcrae82(numberf)
        else:
            raise Exception("Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!")
    
        # Get the names of the parameters to analyse
        names = ["par" + name for name in self.parameter()["name"]]

        # open the file for sensitivity results
        f = open(dbname+".csv", 'w+')
        np.savetxt(f, [names], "%s", delimiter=",")
        
        data = np.array(results)
        
        try:            
            # convert array of void to array of float
            mod_results = np.full((data.shape[0], len(data[0])), np.nan)

            for index in range(len(data)):        
                mod_results[index, :] = list(data[index])
      
            print("data series")
            for index in range(mod_results.shape[1]):    
                x_data = mod_results[:, index]
                sens_data = efast_sensitivity(x_data, numberf, t_runs, t_freq)
                np.savetxt(f, [sens_data], delimiter=",", fmt='%1.5f')
                
        except:
            print("single data")
            x_data = data
            sens_data = efast_sensitivity(x_data, numberf, t_runs, t_freq)
            np.savetxt(f, [sens_data], delimiter=",", fmt='%1.5f')
        
        f.close()

1	# -- coding: utf-8 --
2	"""
3	EFAST algorithm transfered from CRAN R package 'fast' by Dominik Reusser (2015)
4	Copy of most functions and descriptions and implementation in SPOTPY framework by Anna Herzog (2024)
5
6	This Version of FAST differs from the previously impemented Version, as it does not require a precalculation of minimum model runs by the user. It rather uses a hard coded definition of minimum runs based on Cukier or McRae.
7	It therefore requieres considreably less model runs (for 5 Parameters 2245 (fast) vs. 71 runs (efast)).
8
9	For further information on the functions, the origninal R package and the fast algorithm please have a look at the R package description from CRAN, or see the following references:
10
11	Reusser, Dominik E., Wouter Buytaert, and Erwin Zehe.
12	"Temporal dynamics of model parameter sensitivity for computationally expensive models with FAST (Fourier Amplitude Sensitivity Test)." Water Resources Research 47 (2011): W07551.
13
14	Further References:
15	CUKIER, R. I.; LEVINE, H. B. & SHULER, K. E. Non-Linear
16	"Sensitivity Analysis Of Multi-Parameter Model Systems" Journal Of Computational Physics, 1978 , 26 , 1-42
17	CUKIER, R. I.; FORTUIN, C. M.; SHULER, K. E.; PETSCHEK, A. G. & SCHAIBLY, J. H.
18	"Study Of Sensitivity Of Coupled Reaction Systems To Uncertainties In Rate Coefficients .1. Theory" Journal Of Chemical Physics, 1973 , 59 , 3873-3878
19	SCHAIBLY, J. H. & SHULER, K. E.
20	"Study Of Sensitivity Of Coupled Reaction Systems To Uncertainties In Rate Coefficients .2. Applications" Journal Of Chemical Physics, 1973 , 59 , 3879-3888
21	CUKIER, R. I.; SCHAIBLY, J. H. & SHULER, K. E.
22	"Study Of Sensitivity Of Coupled Reaction Systems To Uncertainties In Rate Coefficients .3. Analysis Of Approximations" Journal Of Chemical Physics, 1975 , 63 , 1140-1149
23
24
25	Copyright (c) 2018 by Tobias Houska
26	This file is part of Statistical Parameter Optimization Tool for Python(SPOTPY).
27	:author: Anna Herzog, Tobias Houska, Dominik Reusser
28	"""
29
30	import warnings	9✔
31
32	import numpy as np	9✔
33
34	from ..analyser import efast_sensitivity, get_modelruns	9✔
35	from . import _algorithm	9✔
36
37
38	class efast(_algorithm):	9✔
39	"""
40	EFAST Algorithm for (distributed) parameter Sensitivity after FAST algorithm according to Cukier 1975 or McRae 1982
41	The Fourier Amplitude Sensitivity Test (FAST) is a method to determine global sensitvities of a model on parameter changes
42	within relatively few model runs.
43	"""
44
45	_unaccepted_parameter_types = ()	9✔
46
47	def __init__(self, args, *kwargs):	9✔
48	"""
49	Input
50	----------
51	spot_setup: class
52	model: function
53	Should be callable with a parameter combination of the parameter-function
54	and return an list of simulation results (as long as evaluation list)
55	parameter: function
56	When called, it should return a random parameter combination. Which can
57	be e.g. uniform or Gaussian
58	objectivefunction: function
59	Should return the objectivefunction for a given list of a model simulation and
60	observation.
61	evaluation: function
62	Should return the true values as return by the model.
63
64	dbname: str
65	* Name of the database where parameter, objectivefunction value and simulation results will be saved.
66
67	dbformat: str
68	* ram: fast suited for short sampling time. no file will be created and results are saved in an array.
69	* csv: A csv file will be created, which you can import afterwards.
70
71	parallel: str
72	* seq: Sequentiel sampling (default): Normal iterations on one core of your cpu.
73	* mpi: Message Passing Interface: Parallel computing on cluster pcs (recommended for unix os).
74
75	save_sim: boolean
76	* True: Simulation results will be saved
77	* False: Simulation results will not be saved
78	"""
79
80
81
82	kwargs["algorithm_name"] = "EFast Sampler"	9✔
83	super(efast, self).__init__(args, *kwargs)	9✔
84
85	# hard coded repetition and frequency values from R package Fast
86	d_m_cukier75 = [4,8,6,10,20,22,32,40,38,26,56,62,46,76,96,60,86,126,134,112,92,128,154,196,34,416,106,208,328,198,382,88,348,186,140,170,284,568,302,438,410,248,448,388,596,217,100,488,166]	9✔
87	min_runs_cukier75 = [np.nan, np.nan, 19, 39, 71, 91, 167, 243, 315, 403, 487, 579, 687, 907, 1019, 1223, 1367, 1655, 1919, 2087, 2351, 2771, 3087, 3427, 3555, 4091, 4467, 4795, 5679, 5763, 6507, 7103, 7523, 8351, 9187, 9667, 10211, 10775, 11339, 7467, 12891, 13739, 14743, 15743, 16975, 18275, 18927, 19907, 20759, 21803]	9✔
88	omega_m_cukier75 = [np.nan, np.nan, 1, 5, 11, 1, 17, 23, 19, 25, 41, 31, 23, 87, 67, 73, 58, 143, 149, 99, 119, 237, 267, 283, 151, 385, 157, 215, 449, 163, 337, 253, 375, 441, 673, 773, 875, 873, 587, 849, 623, 637, 891, 943, 1171, 1225, 1335, 1725, 1663, 2019]	9✔
89	d_m_mcrae82 = [4, 8, 6, 10, 20, 22, 32, 40, 38, 26, 56, 62, 46, 76, 96]	9✔
90	min_runs_mcrae82 = [0, 15, 27, 47, 79, 99, 175, 251, 323, 411, 495, 587, 695, 915, 1027]	9✔
91	omega_m_mcrae82 = [0, 3, 1, 5, 11, 1, 17, 23, 19, 25, 41, 31, 23, 87, 67]	9✔
92
93	def freq_cukier(self, m, i = 1, omega_before = -1):	9✔
94
95	"""
96	This function creates an array of independant frequencies accroding to
97	Cukier 1975 for usage in the eFAST method.
98
99	Parameters
100	----------
101	m: int
102	number of parameters (frequencies) needed
103	i: (intern) int
104	internally used recursion counter
105	omega_before: (intern) int
106	internally used previous frequency
107
108	Raises
109	----------
110	Exeption:
111	"Parameter number m to high, not implemented"
112	Execution is stopped, if the number of Parameters is to high (max. Number of Parameters for Cukier = 50)
113
114	Returns
115	----------
116	np.array of float:
117	A 1d-Array of independant frequencies to the order of 4
118
119	Author
120	----------
121	Dominic Reusser
122	spotpy translation: Anna Herzog
123
124	"""
125
NEW 126	if i <= 1:	×
NEW 127	if m >= len(self.omega_m_cukier75):	×
NEW 128	raise Exception("Parameter number m to high, not implemented")	×
129	else:
NEW 130	o = self.omega_m_cukier75[m-1]	×
NEW 131	return np.append(o, self.freq_cukier(m, i+1, o))	×
132	else:
NEW 133	o = omega_before + self.d_m_cukier75[m-i]	×
NEW 134	if i == m:	×
NEW 135	return o	×
136	else:
NEW 137	return np.append(o, self.freq_cukier(m, i+1, o))	×
138
139
140	def freq_mcrae82(self, m, i = 1, omega_before = -1):	9✔
141
142	"""
143	This function creates an array of independant frequencies accroding to
144	McRae 1982 for usage in the eFAST method.
145
146	Parameters
147	----------
148	m: int
149	number of parameters (frequencies) needed
150	i: (intern) int
151	internally used recursion counter
152	omega_before: (intern) int
153	internally used previous frequency
154
155	Raises
156	----------
157	Exeption:
158	"Parameter number m to high, not implemented"
159	Execution is stopped, if the number of Parameters is to high (max. Number of Parameters for McRae = 15)
160
161	Returns
162	----------
163	np.array of float:
164	A 1d-Array of independant frequencies to the order of 4
165
166	Author
167	----------
168	Dominic Reusser
169	spotpy translation: Anna Herzog
170
171	"""
172
NEW 173	if i <= 1:	×
NEW 174	if m >= len(self.omega_m_mcrae82):	×
NEW 175	raise Exception("Parameter number m to high, not implemented")	×
176	else:
NEW 177	o = self.omega_m_mcrae82[m-1]	×
NEW 178	return np.append(o, self.freq_mcrae82(m, i+1, o))	×
179	else:
NEW 180	o = omega_before + self.d_m_mcrae82[m-i]	×
NEW 181	if i == m:	×
NEW 182	return o	×
183	else:
NEW 184	return np.append(o, self.freq_mcrae82(m, i+1, o))	×
185
186
187	def s(self, m, freq = 'cukier'):	9✔
188
189	"""
190	Function that generates a number of equally spaced values between -p1/2
191	and pi/2. The number is determined by the number of runs required of the
192	eFAST method for a number of parameters (min_runs_cukier or min_runs_mcrae)
193
194	Parameters
195	----------
196	m: int
197	number of parameters/ frequencies
198
199	freq: (optional) str
200	indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
201	Default is Cukier
202
203	Raises
204	----------
205	Exeption:
206	"Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!"
207	Execution is stopped if an invalid method for the frequency selection is provided
208
209	Returns
210	----------
211	2D array of float:
212	an array of equally spaced values between -pi/2 and pi/2
213
214	Author
215	----------
216	Dominic Reusser
217	spotpy translation: Anna Herzog
218
219	"""
220
NEW 221	if freq == 'cukier':	×
NEW 222	min_runs = self.min_runs_cukier75	×
NEW 223	elif freq == 'mcrae':	×
NEW 224	min_runs = self.min_runs_mcrae82	×
225	else:
NEW 226	raise Exception("Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!")	×
227
NEW 228	r = min_runs[m-1]	×
NEW 229	r_range = np.array(range(1,r+1))	×
NEW 230	s = np.pi/r * (2*r_range-r-1)/2	×
231
NEW 232	return s	×
233
234
235	def S(self, m, freq = 'cukier'):	9✔
236
237	"""
238	Function to generate an array of values, which provides the base for the parameterdistribution
239	in the FAST method. It is usally not used directly but called from the
240	fast_parameters function.
241
242	Parameters
243	----------
244	m: int
245	number of parameters/frequencies
246
247	freq: (optional) str
248	indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
249	Default is Cukier
250
251	Raises
252	----------
253	Exeption:
254	"Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!"
255	Execution is stopped if an invalid method for the frequency selection is provided
256
257	Returns
258	----------
259	2D array of float:
260	an array with the shape (number of runs, parameters)
261
262	Author
263	----------
264	Dominic Reusser
265	spotpy translation: Anna Herzog
266	"""
267
NEW 268	if freq == 'cukier':	×
NEW 269	omega = self.freq_cukier(m)	×
NEW 270	elif freq == 'mcrae':	×
NEW 271	omega = self.freq_mcrae82(m)	×
272	else:
NEW 273	raise Exception("Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!")	×
274
NEW 275	tab = np.outer(self.s(m, freq), omega)	×
276
NEW 277	toreturn = np.arcsin(np.sin(tab))/np.pi	×
278
279	# naming array dimensions is not possible with numpy but convention would be toreturn.shape () = (runs, parameternames)
280
NEW 281	return toreturn	×
282
283
284	def rerange(self, data, min_goal = 0, max_goal = 1, center = np.nan):	9✔
285
286	"""
287	This function performes a linear transformation of the data, such that
288	afterwards range(data) = (theMin, theMax)
289
290	Parameters
291	----------
292	data: array of flaot
293	an array with the data to transform
294	in this case the parameter distribution generated by function S
295	min_goal: float
296	the new minimum value (lower parameter bound)
297	max_goal: float
298	the new maximum value (upper parameter bound)
299	center: (optional) float
300	indicates which old value should become the new center
301	default: (max_goal+min_goal)/2
302
303	Returns
304	----------
305	2D array of float:
306	array with the transformed data
307
308	Author
309	----------
310	Dominic Reusser
311	spotpy translation: Anna Herzog
312
313	"""
314
NEW 315	min_data = min(data)	×
NEW 316	max_data = max(data)	×
317
NEW 318	if np.isnan(center):	×
NEW 319	max_data = max_data	×
NEW 320	dat = data - min_data	×
NEW 321	dat = dat/(max_data-min_data) * (max_goal-min_goal)	×
NEW 322	dat = dat + min_goal	×
NEW 323	return(dat)	×
324	else:
325	# split linear transformation, the data is split into two segments, one blow center, one above,
326	# each segment undergoes it's own linear transformation
NEW 327	below = data <= center	×
NEW 328	above = data > center	×
NEW 329	new_data = np.copy(data)	×
NEW 330	np.place(new_data, below, self.rerange(np.insert(data[below], 0, center), min_goal = min_goal, max_goal= max_goal/2)[1:])	×
NEW 331	np.place(new_data, above, self.rerange(np.insert(data[above], 0, center), min_goal = max_goal/2, max_goal= max_goal)[1:])	×
NEW 332	return(new_data)	×
333
334
335	def fast_parameters(self, minimum, maximum, freq = 'cukier', logscale = np.nan, names = np.nan):	9✔
336
337	"""
338	Function for the creation of a FAST Parameter set based on thier range
339
340	Parameters
341	----------
342	minimum: array of float
343	array containing the lower parameter boundries
344	maximum: array of float
345	array containing the upper parameter boundries
346	names: (optional) str
347	array containing the parameter names
348	freq: (optional) str
349	indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
350	Default is Cukier
351	logscale: (optional) bool
352	array containing bool values indicating weather a parameter is varied
353	on a logarithmic scale. In that case minimum and maximum are exponents
354
355	Raises
356	----------
357	Exeption:
358	"Expecting minimum and maximum of same size"
359	Execution is stopped if the number of minimum, maximum or logscale values differ
360
361	Returns
362	----------
363	value: array of float
364	an array with the shape (number of runs, parameters) containing the
365	required parameter sets
366
367	Author
368	----------
369	Dominic Reusser
370	spotpy translation: Anna Herzog
371
372	"""
373
NEW 374	if np.isnan(logscale):	×
NEW 375	logscale = np.full(minimum.shape, False)	×
NEW 376	if np.isnan(names):	×
NEW 377	names = ["P"+str(i) for i in range(minimum.shape[0])]	×
NEW 378	names = np.array(names)	×
379
NEW 380	n = len(minimum)	×
381
NEW 382	if (n != len(maximum)):	×
NEW 383	raise Exception("Expecting minimum and maximum of same size")	×
NEW 384	elif(n != len(names)):	×
NEW 385	raise Exception("Expecting minimum and names of same size")	×
NEW 386	elif(n != len(logscale)):	×
NEW 387	raise Exception("Expecting minimum and logscale of same size")	×
388
NEW 389	toreturn = self.S(m=n, freq = freq) #par_names = names,	×
390
NEW 391	for i in range(0,n):	×
NEW 392	toreturn[:,i] = self.rerange(toreturn[:,i], minimum[i], maximum[i])	×
NEW 393	if logscale[i]:	×
NEW 394	toreturn[:,i] = 10**toreturn[:,i]	×
395
NEW 396	return toreturn	×
397
398
399	def sample(self, repetitions, freq = 'cukier', logscale = np.nan):	9✔
400	"""
401	Samples from the eFAST algorithm.
402
403	Input
404	----------
405	repetitions: int
406	Maximum number of runs.
407	freq: (optional) str
408	indicates weather to use the frequencies after 'cukier' or McRae 'mcrae'
409	Default is Cukier
410	logscale: (optional) bool
411	array containing bool values indicating weather a parameter is varied
412	on a logarithmic scale. In that case minimum and maximum are exponents
413
414	Raises
415	----------
416	Warning:
417	"specified number of repetitions is smaller than minimum required for eFAST analysis!
418	Simultions will be perfomed but eFAST analysis might not be possible"
419
420	If the specified number of repetitions is smaller than required by the algorithm, the
421	algortihm will still execute, but the usage of spotpy.analyser.efast.calc_sensitivity()
422	might not be possible. Rather define a large number of reps as only required minimum runs
423	will be performed anyway.
424
425	Warning:
426	"Number of specified repetitions equals or exeeds number of minimum required repetitions for eFAST analysis
427	programm will stop after required runs."
428
429	Returns
430	----------
431	database:
432	spotpy data base in chosen format with objective function values, parameter values and simulation results
433
434	Author
435	----------
436	Tobias Houska, Anna Herzog
437	"""
438
NEW 439	print("generating eFAST Parameters")	×
440	# Get the names of the parameters to analyse
NEW 441	names = self.parameter()["name"]	×
442	# Get the minimum and maximum value for each parameter from the
443	# distribution
NEW 444	parmin, parmax = self.parameter()["minbound"], self.parameter()["maxbound"]	×
NEW 445	self.numberf = len(parmin)	×
NEW 446	min_runs = self.min_runs_cukier75[self.numberf-1]	×
447
NEW 448	if min_runs > repetitions:	×
NEW 449	warnings.warn("specified number of repetitions is smaller than minimum required for eFAST analysis!\n"+	×
450	"Simultions will be perfomed but eFAST analysis might not be possible")
451	else:
NEW 452	repetitions = min_runs	×
NEW 453	print("Number of specified repetitions equals or exeeds number of minimum required repetitions for eFAST analysis\n"+	×
454	"programm will stop after required runs.")
455
NEW 456	self.set_repetiton(repetitions)	×
457
458	# Generate Matrix with eFAST parameter sets
NEW 459	N = self.fast_parameters(parmin, parmax, freq, logscale)	×
460
NEW 461	print("Starting the eFast algorithm with {} repetitions...".format(repetitions))	×
462
463	# A generator that produces parametersets if called
NEW 464	param_generator = (	×
465	(rep, N[rep,:]) for rep in range(int(repetitions))
466	)
NEW 467	for rep, randompar, simulations in self.repeat(param_generator):	×
468	# A function that calculates the fitness of the run and the manages the database
NEW 469	self.postprocessing(rep, randompar, simulations)	×
470
NEW 471	if self.breakpoint == "write" or self.breakpoint == "readandwrite":	×
NEW 472	if rep >= lastbackup + self.backup_every_rep:	×
NEW 473	work = (rep, N[rep,:])	×
NEW 474	self.write_breakdata(self.dbname, work)	×
NEW 475	lastbackup = rep	×
NEW 476	self.final_call()	×
477
478	def calc_sensitivity(self, results, dbname, freq = 'cukier'):	9✔
479	"""
480	Function to calcultae the sensitivity for a data series (eg. a time series)
481	based on the eFAST algorithm.
482
483	Parameters
484	----------
485	results: np.array of void
486	spopty restults array as loaded with spotpy.analyser.load_csv_results()
487	dbname: (optional) str
488	name of file to store sensitivity values
489	freq: (optional) str
490	indicates weather to use the frequencies after Cukier 'cukier' or McRae 'mcrae'
491	Default is Cukier
492
493	Raises
494	----------
495	Exeption:
496	"Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!"
497	Execution is stopped if an invalid method for the frequency selection is provided
498
499	Returns
500	----------
501	database:
502	database containing the temporal parameter sensitivities with seperate column
503	for each parameter and row for eg. time step
504
505	Author
506	----------
507	Dominic Reusser
508	spotpy translation: Anna Herzog
509
510	"""
NEW 511	numberf = len(self.parameter()["minbound"])	×
512
513	# sens_data = np.full((len(results[0]), numberf), np.nan)
514
515	# idendtify frequency
NEW 516	if freq == 'cukier':	×
NEW 517	t_runs = self.min_runs_cukier75[numberf-1]	×
NEW 518	t_freq = self.freq_cukier(numberf)	×
NEW 519	elif freq == 'mcrae':	×
NEW 520	t_runs = self.min_runs_mcrae82[numberf-1]	×
NEW 521	t_freq = self.freq_mcrae82(numberf)	×
522	else:
NEW 523	raise Exception("Missing option for Frequency selection for Parameter definition, choose cukier or mcrae!")	×
524
525	# Get the names of the parameters to analyse
NEW 526	names = ["par" + name for name in self.parameter()["name"]]	×
527
528	# open the file for sensitivity results
NEW 529	f = open(dbname+".csv", 'w+')	×
NEW 530	np.savetxt(f, [names], "%s", delimiter=",")	×
531
NEW 532	data = np.array(results)	×
533
NEW 534	try:	×
535	# convert array of void to array of float
NEW 536	mod_results = np.full((data.shape[0], len(data[0])), np.nan)	×
537
NEW 538	for index in range(len(data)):	×
NEW 539	mod_results[index, :] = list(data[index])	×
540
NEW 541	print("data series")	×
NEW 542	for index in range(mod_results.shape[1]):	×
NEW 543	x_data = mod_results[:, index]	×
NEW 544	sens_data = efast_sensitivity(x_data, numberf, t_runs, t_freq)	×
NEW 545	np.savetxt(f, [sens_data], delimiter=",", fmt='%1.5f')	×
546
NEW 547	except:	×
NEW 548	print("single data")	×
NEW 549	x_data = data	×
NEW 550	sens_data = efast_sensitivity(x_data, numberf, t_runs, t_freq)	×
NEW 551	np.savetxt(f, [sens_data], delimiter=",", fmt='%1.5f')	×
552
NEW 553	f.close()	×

thouska / spotpy / 13288720301

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous