Skip to content
Snippets Groups Projects
Commit ab243150 authored by Philipp Frank's avatar Philipp Frank
Browse files

Merge branch 'misc_update' into 'main' 

Misc update

See merge request !7
parents a28324c5 60e5ddc5
No related branches found
No related tags found
1 merge request!7Misc update
Pipeline #149176 passed
Stage: test
Show whitespace changes
Inline Side-by-side
.gitlab-ci.yml
+ 1
− 1
View file @ ab243150
......@@ -6,5 +6,5 @@ run:
- pip install matplotlib nifty8 ift-resolve
- jupyter nbconvert --execute --ExecutePreprocessor.timeout=None --to html demo_CorrelatedFields.ipynb
- jupyter nbconvert --execute --ExecutePreprocessor.timeout=None --to html demo_radio.ipynb
- jupyter nbconvert --execute --ExecutePreprocessor.timeout=None --to html demo_poisson.ipynb
- jupyter nbconvert --execute --ExecutePreprocessor.timeout=None --to html demo_wiener_filter.ipynb
demo_poisson.ipynb
+ 30
− 10
View file @ ab243150
%% Cell type:markdown id: tags:
Nifty tutorial for Poisson count data
=====================================
%% Cell type:markdown id: tags:
Setup
-----
%% Cell type:code id: tags:
``` python
import numpy as np
import matplotlib.pyplot as plt
import nifty8 as ift
from utils import plot_2D, load_psf, geovi_sampling, plot_posterior
ift.random.push_sseq_from_seed(42)
evidences = np.zeros(3)
```
%% Cell type:code id: tags:
``` python
# Load data and visualize
data = np.load('data/poisson.npz')
print(data['data'].shape)
plot_2D(data['data'], 'Data')
```
%% Cell type:code id: tags:
``` python
position_space = ift.RGSpace([128, 128])
# Homogeneous poisson process
projection = ift.VdotOperator(ift.full(position_space, 1.)).adjoint
model1 = ift.FieldAdapter(projection.domain, 'hom')
model1 = ift.exp(5. * model1)
model1 = projection @ model1
print(model1)
```
%% Cell type:code id: tags:
``` python
# Set up likelihood & PSF
d = ift.makeField(position_space, data['data'])
likelihood = ift.PoissonianEnergy(d)
PSF_op, psf = load_psf(position_space)
plot_2D(psf, 'PSF')
likelihood = likelihood @ PSF_op
```
%% Cell type:code id: tags:
``` python
# Inference model 1
samples, evidence = geovi_sampling(likelihood @ model1)
plot_2D(samples.average(model1).val, 'model1 posterior mean')
print(evidence)
evidences[0] = evidence
```
%% Cell type:code id: tags:
``` python
# Independent poisson process
model2 = ift.InverseGammaOperator(position_space, 2., 3.).ducktape('independent')
# Inference model 2
samples, evidence = geovi_sampling(likelihood @ model2)
plot_2D(samples.average(model2).val, 'model2 posterior mean')
print(evidence)
evidences[1] = evidence
```
%% Cell type:code id: tags:
``` python
# Compare evidence
print(evidences[:2])
```
%% Cell type:code id: tags:
``` python
# Diffuse poisson process
args = {
'offset_mean': .5,
'offset_std': (1., 1E-5),
# Amplitude of field fluctuations
'fluctuations': (1.5, 0.5), # 1.0, 1e-2
# Exponent of power law power spectrum component
'loglogavgslope': (-4., 1), # -6.0, 1
# Amplitude of integrated Wiener process power spectrum component
'flexibility': (1., 0.2), # 2.0, 1.0
# How ragged the integrated Wiener process component is
'asperity': (0.1, 0.01), # 0.1, 0.5
# Name of the input keys
'prefix' : 'diffuse'
}
correlated_field = ift.SimpleCorrelatedField(position_space, **args)
pspec = correlated_field.power_spectrum
diffuse = correlated_field.exp()
model3 = diffuse + model2
```
%% Cell type:code id: tags:
``` python
# Prior samples
pl = ift.Plot()
for _ in range(9):
pl.add(model3(ift.from_random(model3.domain)))
pl.output()
```
%% Cell type:code id: tags:
``` python
# Inference model 3
samples, evidence = geovi_sampling(likelihood @ model3)
plot_2D(samples.average(model3).val, 'model3 posterior mean')
print(evidence)
evidences[2] = evidence
```
%% Cell type:code id: tags:
``` python
# Compare evidence
print(evidences)
print(evidences - evidences[-1])
```
%% Cell type:markdown id: tags:
Posterior visualization
-----------------------
%% Cell type:code id: tags:
``` python
plot_posterior(samples, data, model3, diffuse, model2, pspec)
```
......
demo_poisson_solution.ipynb deleted 100644 → 0
+ 0
− 230
View file @ a28324c5
%% Cell type:markdown id: tags:
Nifty tutorial for Poisson count data
=====================================
%% Cell type:markdown id: tags:
Setup
-----
%% Cell type:code id: tags:
``` python
import numpy as np
import matplotlib.pyplot as plt
import nifty8 as ift
from utils import plot_2D, load_psf, geovi_sampling, plot_posterior
ift.random.push_sseq_from_seed(42)
evidences = np.zeros(3)
```
%% Cell type:code id: tags:
``` python
# Load data and visualize
data = np.load('data/poisson.npz')
print(data['data'].shape)
plot_2D(data['data'], 'Data')
```
%% Cell type:code id: tags:
``` python
position_space = ift.RGSpace([128, 128])
# Homogeneous poisson process
projection = ift.VdotOperator(ift.full(position_space, 1.)).adjoint
model1 = ift.FieldAdapter(projection.domain, 'hom')
model1 = ift.exp(5. * model1)
model1 = projection @ model1
print(model1)
```
%% Cell type:code id: tags:
``` python
# Set up likelihood & PSF
d = ift.makeField(position_space, data['data'])
likelihood = ift.PoissonianEnergy(d)
PSF_op, psf = load_psf(position_space)
plot_2D(psf, 'PSF')
likelihood = likelihood @ PSF_op
```
%% Cell type:code id: tags:
``` python
# Inference model 1
samples, evidence = geovi_sampling(likelihood @ model1)
plot_2D(samples.average(model1).val, 'model1 posterior mean')
print(evidence)
evidences[0] = evidence
```
%% Cell type:code id: tags:
``` python
# Independent poisson process
model2 = ift.InverseGammaOperator(position_space, 2., 3.).ducktape('independent')
# Inference model 2
samples, evidence = geovi_sampling(likelihood @ model2)
plot_2D(samples.average(model2).val, 'model2 posterior mean')
print(evidence)
evidences[1] = evidence
```
%% Cell type:code id: tags:
``` python
# Compare evidence
print(evidences[:2])
```
%% Cell type:code id: tags:
``` python
# Diffuse poisson process
args = {
'offset_mean': .5,
'offset_std': (1., 1E-5),
# Amplitude of field fluctuations
'fluctuations': (1.5, 0.5), # 1.0, 1e-2
# Exponent of power law power spectrum component
'loglogavgslope': (-4., 1), # -6.0, 1
# Amplitude of integrated Wiener process power spectrum component
'flexibility': (1., 0.2), # 2.0, 1.0
# How ragged the integrated Wiener process component is
'asperity': (0.1, 0.01), # 0.1, 0.5
# Name of the input keys
'prefix' : 'diffuse'
}
correlated_field = ift.SimpleCorrelatedField(position_space, **args)
pspec = correlated_field.power_spectrum
diffuse = correlated_field.exp()
model3 = diffuse + model2
```
%% Cell type:code id: tags:
``` python
# Prior samples
pl = ift.Plot()
for _ in range(9):
pl.add(model3(ift.from_random(model3.domain)))
pl.output()
```
%% Cell type:code id: tags:
``` python
# Inference model 3
samples, evidence = geovi_sampling(likelihood @ model3)
plot_2D(samples.average(model3).val, 'model3 posterior mean')
print(evidence)
evidences[2] = evidence
```
%% Cell type:code id: tags:
``` python
# Compare evidence
print(evidences)
print(evidences - evidences[-1])
```
%% Cell type:markdown id: tags:
Posterior visualization
-----------------------
%% Cell type:code id: tags:
``` python
plot_posterior(samples, data, model3, diffuse, model2, pspec)
```
variational_inference_visualized.py
+ 1
− 1
View file @ ab243150
......@@ -36,7 +36,7 @@ def main():
dom = ift.UnstructuredDomain(1)
n_samples = 20
show_geo = False
show_geo = True
scale = 10.
def transformation(x,y):
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment