Merge pull request #7 from tudo-physik-e4/dev

merge Dev

.gitignore

@@ -7,3 +7,4 @@
 .vscode
 
 Manifest.toml
+*.pdf

Project.toml

@@ -1,23 +1,28 @@
 name = "EFTfitter"
 uuid = "2a2a6a96-a4ec-477f-941b-476720990667"
 authors = ["Cornelius Grunwald"]
-version = "0.1.2"
+version = "0.2"
 
 [deps]
 BAT = "c0cd4b16-88b7-57fa-983b-ab80aecada7e"
+BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d"
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 NamedTupleTools = "d9ec5142-1e00-5aa0-9d6a-321866360f50"
 Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a"
+Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
 Requires = "ae029012-a4dd-5104-9daa-d747884805df"
 Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"
+SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
+TypedTables = "9d95f2ec-7b3d-5a63-8d20-e2491e220bb9"
 ValueShapes = "136a8f8c-c49b-4edb-8b98-f3d64d48be8f"
 
 [compat]
-BAT = "2, 3"
+BAT = "3"
 Distributions = "0.22, 0.23, 0.24, 0.25"
 NamedTupleTools = "0.13"
 Parameters = "0.12"
@@ -26,7 +31,7 @@ Requires = "0.5, 1"
 Setfield = "0.7"
 StatsBase = "0.32, 0.33"
 ValueShapes = "0.7, 0.8, 0.9, 0.10"
-julia = "1.3"
+julia = "1.8, 1.9"
 
 [extras]
 IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953"

README.md

@@ -10,11 +10,18 @@
 
 New implementation of the [EFTfitter](https://github.com/tudo-physik-e4/EFTfitterRelease) in the [Julia languange](https://julialang.org/).  
 
-Tool for constraining the parameters of physics models using Bayesian inference by combining measurements of (different) observables.
-Particularly suited for EFT (effective field theory) interpretations. 
+EFTfitter is a tool for constraining the parameters of physics models using Bayesian inference by combining measurements of (different) observables.
+It is particularly suited for EFT (effective field theory) interpretations, but is not limited to these use cases.
 
 Work-in-progress, interfaces and functionalities might be subject to changes.
 
+## News
+Since version 0.2 of the package, the following new features are available:
+- Model uncertainties: The functions giving the predictions for the observable values can now also return a parameter-dependent value quantifying the uncertainty on the prediction. These uncertainties are currently treated as uncorrelated and are added to the total covariance matrix.
+- The data type of the (inverse) covariance matrix can now be set by the user. This can allow to increase the performance by speeding up the vector-matrix multiplication, e.g. in the case of sparse covariance matrices.
+- `MeasurementDistribution` is now called `BinnedMeasurement`
+
+
 ## Installation
 The EFTfitter.jl package can be installed using:
 ```julia
@@ -27,7 +34,7 @@ Please see the [installation guide](https://tudo-physik-e4.github.io/EFTfitter.j
 
 ## Documentation & Tutorials
 Please see the [documentation](https://tudo-physik-e4.github.io/EFTfitter.jl/dev/) for tutorials and information on how to use EFTfitter.jl.  
-Executable versions of the tutorials and examples can also be found [here](https://github.com/tudo-physik-e4/EFTfitter.jl/tree/main/examples/tutorials).
+Executable versions of the tutorials and examples can also be found [here](https://github.com/tudo-physik-e4/EFTfitter.jl/tree/main/examples).
 
 You can also try running the tutorials right now: [![badge](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/tudo-physik-e4/EFTfitter.jl/binder?urlpath=git-pull%3Frepo%3Dhttps%253A%252F%252Fgithub.com%252Ftudo-physik-e4%252FEFTfitter.jl%26urlpath%3Dtree%252FEFTfitter.jl%252Fexamples%252Fnotebooks%252F%26branch%3Dmain)
 

docs/src/BLUE.md

@@ -10,45 +10,45 @@ by L. Lyons, D. Gibaut and P. Clifford (https://www.sciencedirect.com/science/ar
 All numbers are taken from the example on charm particle lifetime experiments in section 5.
 A factor of 10^13 is applied for convenience.
 
-```julia
+````julia
 using EFTfitter
 using BAT
 using IntervalSets
 using Statistics
 using StatsBase
 using LinearAlgebra
 using Plots
-```
+````
 
 We need one parameter for the best estimator and choose
 a uniform distribution in the range 8 to 14 as prior:
 
-```julia
+````julia
 parameters = BAT.NamedTupleDist(
     τ = 8..14,
 )
-```
+````
 
 When combining multiple measurements of the same observable,
 only a function returning the combination parameter is needed:
 
-```julia
+````julia
 estimator(params) = params.τ
-```
+````
 
 In Eq. (17') of the reference paper the following covariance matrix is given:
 
-```julia
+````julia
 covariance = [2.74 1.15 0.86 1.31;
               1.15 1.67 0.82 1.32;
               0.86 0.82 2.12 1.05;
               1.31 1.32 1.05 2.93]
-```
+````
 
 For using this in EFTfitter.jl, we first need to convert the covariance matrix
 into a correlation matrix and the corresponding uncertainty values:
 
-```julia
+````julia
 corr, unc = EFTfitter.cov_to_cor(covariance)
 
 measurements = (
@@ -61,52 +61,52 @@ measurements = (
 correlations = (
     stat = Correlation(corr),
 )
-```
+````
 
 construct an `EFTfitterModel`:
 
-```julia
+````julia
 model = EFTfitterModel(parameters, measurements, correlations)
 posterior = PosteriorMeasure(model);
-```
+````
 
 sample the posterior with BAT.jl:
 
-```julia
+````julia
 algorithm = MCMCSampling(mcalg =MetropolisHastings(), nsteps = 10^6, nchains = 4)
 samples = bat_sample(posterior, algorithm).result
-```
+````
 
 plot the posterior distribution for the combination parameter τ:
 
-```julia
+````julia
 plot(samples, :τ, mean=true)
-```
+````
 
 ![blue plots](plots/plot_blue.png)
 
 print numerical results of combination:
 
-```julia
+````julia
 println("Mode: $(mode(samples).τ)")
 println("Mean: $(mean(samples).τ) ± $(std(samples).τ)")
 ```
 Mode: 11.15985
 Mean: 11.15471 ± 0.80180
 ```
-```
+````
 
 ### Comparison with BLUE method
 
-```julia
+````julia
 blue = BLUE(model)
 println("BLUE: $(blue.value) ± $(blue.unc)")
 println("BLUE weights: $(blue.weights)")
 ```
 BLUE: 11.15983 ± 1.28604
 BLUE weights: [0.145, 0.470, 0.347, 0.038]
 ```
-```
+````
 
 ---
 

docs/src/advanced_tutorial.md

@@ -6,20 +6,20 @@ Pages = ["advanced_tutorial.md"]
 Depth = 3
 ```
 
-## Vector of functions for a MeasurementDistribution
-When using distributions of measurements, a vector of functions with the predictions
-for the observable needs to be passed containing a function for each of the bins which
-have only the model parameters as their argument. Defining a separate function for each
+## Vector of functions for a BinnedMeasurement
+When using binned measurements, a vector of functions giving the predictions
+for the observable needs to be passed. It contains a function for each of bin and
+has only the model parameters as its argument. Defining a separate function for each
 bin can, however, become tedious for a large number of bins, especially since typically
 the bins of a distribution have a similar functional dependence on the model parameters
 and only differ in some coefficients. In such cases, it is possible to use Julia's
-[metaprogramming](https://docs.julialang.org/en/v1/manual/metaprogramming/) features to
-create the vector of functions. The distribution in our [basic tutorial](https://tudo-physik-e4.github.io/EFTfitter.jl/dev/tutorial/)
+anonymous functions to quickly create the vector of functions.
+The distribution in our [basic tutorial](https://tudo-physik-e4.github.io/EFTfitter.jl/dev/tutorial/)
 has been defined by implementing three functions that all call the same function `myfunc`
 but with different values for the coefficients
 The same result can also be achieved like this:
 
-```julia
+````julia
 function get_coeffs(i) # return the coefficients for bin i
     coeffs = [[2.2, 5.5, 6.6], [2.2, 5.5, 6.6], [2.2, 5.5, 6.6]]
     return coeffs[i]
@@ -29,21 +29,13 @@ function my_dist_func(params, i)
     coeffs = get_coeffs(i)
     return coeffs[1] * params.C1 + coeffs[2] * params.C1 * params.C2+ coeffs[3] * params.C2
 end
-```
+````
 
-create an array of Functions with names `diff_xsec_binX`:
-
-```julia
-diff_xsec=Function[]
-for i in 1:3
-    @eval begin
-        function $(Symbol("diff_xsec_bin$i"))(params)
-            return my_dist_func(params, $i)
-        end
-        push!(diff_xsec, $(Symbol("diff_xsec_bin$i")))
-    end
-end
-```
+create an array of anonymous functions
+
+````julia
+diff_xsec = Function[x -> my_dist_func(x, i) for i in 1:3]
+````
 
 ## Using covariance matrices
 Information about the uncertainties of measurements need to be provided to EFTfitter.jl
@@ -53,17 +45,17 @@ to correlation matrices and uncertainty values before. The function
 [`cov_to_cor`](https://tudo-physik-e4.github.io/EFTfitter.jl/dev/api/#EFTfitter.cov_to_cor-Tuple{Array{var%22#s58%22,2}%20where%20var%22#s58%22%3C:Real})
 can be used for this:
 
-```julia
+````julia
 cov_syst = [3.24   0.81   0.378  0.324  0.468;
             0.81   0.81   0.126  0.162  0.234;
             0.378  0.126  0.49   0.126  0.182;
             0.324  0.162  0.126  0.81   0.234;
             0.468  0.234  0.182  0.234  1.69]
 
 cor_syst, unc_syst = cov_to_cor(cov_syst)
-```
+````
 
-```julia
+````julia
 measurements = (
 
     Meas1 = Measurement(xsec1, 21.6,
@@ -73,23 +65,23 @@ measurements = (
             uncertainties = (stat=0.6, syst=unc_syst[2], another_unc=1.1), active=true),
 
 
-    MeasDist = MeasurementDistribution(diff_xsec, [1.9, 2.93, 4.4],
+    MeasDist = BinnedMeasurement(diff_xsec, [1.9, 2.93, 4.4],
                 uncertainties = (stat = [0.7, 1.1, 1.2], syst= unc_syst[3:5], another_unc = [1.0, 1.2, 1.9]),
                 active=[true, false, true]),
 )
-```
+````
 
 ## Nuisance Correlations
 When performing an analysis with unknown correlation coefficients, it is possible
 to treat them as nuisance parameters in the fit.
 For this, we define a further `NamedTuple` consisting of `NuisanceCorrelation` objects:
 
-```julia
+````julia
 nuisance_correlations = (
     ρ1  = NuisanceCorrelation(:syst, :Meas1, :Meas2, -1..1),
     ρ2  = NuisanceCorrelation(:syst, :MeasDist_bin1, :MeasDist_bin3, truncated(Normal(0.5, 0.1), -1, 1)),
 )
-```
+````
 
 In the `NuisanceCorrelation` object we specify the name of the uncertainty type, the
 names of the two measurements we want to correlate using the nuisance correlations and
@@ -102,12 +94,12 @@ the normal distribution accordingly.
 
 We need to modify the definition of the `EFTfitterModel` by also passing the `nuisance_correlations`:
 
-```julia
-model = EFTfitterModel(parameters, measurements, correlations, nuisance_correlations)
+````julia
+model = EFTfitterModel(parameters, measurements, correlations, nuisances = nuisance_correlations)
 
 
 savefig(p, "plot.pdf")
-```
+````
 
 ![plot with nuisances](plots/plot_nuisance.png)
 !!! warning
@@ -129,28 +121,28 @@ smallest interval containing 90% of the posterior probability when deactivating
 For models with more than one parameter, the sum of the relative increases of all
 one-dimensional smallest intervals is used, i.e. `SumOfSmallestIntervals(p=0.9, bins=200)`.
 
-```julia
+````julia
 measurement_ranking = EFTfitter.rank_measurements(model)
-```
+````
 
 The sampling algorithm to be used can be passed with the keyword `sampling_algorithm`.
 By default, `BAT.MCMCSampling()` is used, i.e. Metropolis-Hastings with 4 chains and 100000 steps.
 
-```julia
+````julia
 plot(measurement_ranking, title = "Ranking of measurements")
-```
+````
 
 ![measurement ranking](plots/meas_ranks.png)
 
 For ranking the uncertainty types, the relative decrease is used.
 
-```julia
+````julia
 uncertainty_ranking = EFTfitter.rank_uncertainties(model,
     criterion = SumOfSmallestIntervals(p=0.9, bins=200),
     sampling_algorithm = SobolSampler(nsamples = 10^5), order = :values)
 
 plot(uncertainty_ranking, title = "Ranking of uncertainty types")
-```
+````
 
 ![measurement ranking](plots/unc_ranks.png)
 Please see the [ranking documentation](https://tudo-physik-e4.github.io/EFTfitter.jl/dev/api/#EFTfitter.rank_measurements) for further ranking criteria and keyword arguments.