[A Matlab script for the following example is avalable at examples/case_tomography/sippi_AM13_metropolis_modeling_error.m.]
The use of any of the forward models defined above, will be approximation to solving the perfect forward problem. This leads to a 'modeling' error as demonstrated by [HCM14]. If one has access to an optimal (but perhaps computational inefficient) forward model, and a faster (less accurate) forward model, then a Gaussian model of the modeling error caused by using the approximate, as opposed to the optima, forward model can be estimated using sippi_compute_modelization_forward_error. sippi_compute_modelization_forward_error.
SIPPI allows accounting for such modeling error through the dt
and Ct
fields for the data structure.
The setup of the data, prior and forward structures is identical to the one described in the previous example.
%% Load the travel time data set from ARRENAES clear all;close all D=load('AM13_data.mat'); options.txt='AM13'; %% SETUP DATA id=1; data{id}.d_obs=D.d_obs; data{id}.d_std=D.d_std; data{id}.Ct=D.Ct+1; % Covariance describing modeling error %% SETUP PRIOR im=1; prior{im}.type='FFTMA'; prior{im}.name='Velocity (m/ns)'; prior{im}.m0=0.145; prior{im}.Va='.0003 Sph(6)'; dx=0.15; prior{im}.x=[-1:dx:6]; prior{im}.y=[0:dx:13]; prior{im}.cax=[.1 .18]; % SETUP THE FORWARD MODEL USED IN INVERSION forward.forward_function='sippi_forward_traveltime'; forward.sources=D.S; forward.receivers=D.R; forward.type='fat';forward.linear=1;forward.freq=0.1;
In order to compute the modeling with respect to using the 'Born' type forward model, one can define a new forward structure, here forward_full
, and estimate a Gaussian model for the modeling error using
% SETUP THE 'OPTIMAL' FORWARD MODEL forward_full.forward_function='sippi_forward_traveltime'; forward_full.sources=D.S; forward_full.receivers=D.R; forward_full.type='Born';forward_full.linear=1;forward_full.freq=0.1; % COMPUTE MODELING ERROR DUE TO USE OF forward AS OPPOSED TO forward_full N=100; [Ct,dt,dd]=sippi_compute_modelization_forward_error(forward_full,forward,prior,data,N); % ASSIGN MODELING ERROR TO DATA data{1}.dt=dt{1 }; data{1}.Ct=Ct{1};
Sampling of the posterior can proceed exactly as for the previous example, using
options.mcmc.nite=500000; % optional, default:nite=30000 options.mcmc.i_sample=500; % optional, default:i_sample=500; options.mcmc.i_plot=1000; % optional, default:i_plot=50; options=sippi_metropolis(data,prior,forward,options); % plot posterior statistics sippi_plot_posterior(options.txt);
This site is hosted by sourceforge.net |