IES phi factors

[cnicol-gwlogic]

Hi there,

When using "ies_phi_factor_file" in pestpp-ies 5.2.4 (and ies_phi_factors_by_real=false), I see the following reported to the console:

Note the two bottom obs groups:
...file tag 'WEL_FLUX' adjusted mean phi (factor): 6.83766e+07 (0.05)
...file tag 'WT_DEPTH_PENALTY' adjusted mean phi (factor): 6.83766e+07 (0.05)

Yet in the phi summary below that phi factor / adjustment reporting, these two are reported as having mean phi values of ~1e+10, or 58% and 18% respectively - not the 6e+7 / 5% of total phi suggested by the factoring process above. This also affects all other obs groups - these two are just an example.

Is this right? Or have I messed something up? I'm thinking the obs noise (which I am applying) could be a factor here - although these two obs groups are penalties with zero noise applied.

Update: I note this in the log file: "using realization BASE residuals for reweighting (if available" - I thought the reweighting was done via the mean across all realisations - no? This could be the cause of the point I identify above, for sure. Can I change that? Maybe I need to turn off "center on base".

Another question for you if I may:
With these numbers being so large, should I scale weights universally down to get to more reasonable numbers? Could the size become an issue numerically?

Thanks in advance.

[cnicol-gwlogic]

Ok, think I've answered my own question - the issue is that I had left "center on base" = true, which I should have turned off anyway for my multimodal update case, per @jtwhite79's suggestion here...if that option is true, then internal reweighting is done on the basis of the base realisation. For multimodal upgrades, this is not relevant, so should not be used for internal reweighting (I would have thought anyway).

It's also worth noting for the record that when "center on base" option was true, and multimodal updgrades are used, I was seeing thread errors on the multimodal upgrade process. When I turn "center on base" == false, these no longer occur. Eg errors (only dumped to console, not logged):

--- starting solve for iteration: 1 ---
...current lambda: 100000
...starting calcs for glm factor 10000
...multimodal solve for inflation factor 10000
Error condition prevents further execution:
thread 6raised an exception: bad allocation thread 7raised an exception: bad allocation thread 8raised an exception: bad allocation thread 9raised an exception: bad allocation thread 10raised an exception: bad allocation thread 11raised an exception: bad allocation thread 12raised an exception: bad allocation thread 13raised an exception: bad allocation thread 14raised an exception: bad allocation thread 15raised an exception: bad allocation thread 16raised an exception: bad allocation thread 17raised an exception: bad allocation thread 18raised an exception: bad allocation thread 19raised an exception: bad allocation thread 20raised an exception: bad allocation thread 21raised an exception: bad allocation thread 22raised an exception: bad allocation thread 23raised an exception: bad allocation thread 24raised an exception: bad allocation thread 25raised an exception: bad allocation thread 26raised an exception: bad allocation thread 27raised an exception: bad allocation thread 28raised an exception: bad allocation thread 29raised an exception: bad allocation thread 30raised an exception: bad allocation thread 31raised an exception: bad allocation

This was what my log file looked like:
22:33:35,2.02e-05,calc initial phi
22:33:36,0.558984,initial phi summary
22:33:36,0.0398565,using full (MAP) update solution
22:33:36,0.0262332,current lambda: 100000
22:33:36,9.23e-05,initialization complete
22:33:36,0.0014942,starting solve for iteration: 1
22:33:36,1.07e-05,starting solve for iteration: 1
22:33:36,6.36e-05,current lambda: 100000
22:33:36,1e-05,reordering variables in pe
22:33:36,0.244311,starting REDSVD
22:33:38,2.54557,triming REDSVD components to 179elements
22:33:39,0.0339867,done REDSVD
22:33:39,0.246896,preparing EnsembleSolver
23:08:03,2064.2,getting phi vectors for all weights
23:08:05,1.58481,starting multithreaded MM neighbor calcs
23:11:34,209.282,starting calcs for glm factor 10000
23:11:34,0.0002514,see .log file for more details
23:11:34,0.0565659,multimodal solve for inflation factor 10000
23:11:34,0.51428,launching threads
23:11:40,5.44262,waiting to join threads
23:15:57,257.082,upgrade thread progress: 100 of 180 parts done
(and that was the end of that file).

Cheers.

[jtwhite79]

Thanks Chris - I'll have a look at that combination of inputs...

[jtwhite79]

The phi factors stuff by default works on the base realization residuals if you passed "center_on_base" or on the mean prior residuals, so essentially it replicates the workflow we are all used to WRT balancing the objective function.

If you want to explore the effect of differing weighting strategies, you have to use a weight ensemble in combination with the multimodal upgrade - this causes each realization to use a unique weight vector (which also effects which realizations that are found in its neighborhood). This can drive each realization to very different minima...

There is also an even crazier option of supplying the phi factors file as a csv file where the index is the same realization names as the noise ensemble and the columns of the "tags" for the different obs (this requires a weight ensemble also) Then each realization should get rebalanced by its residuals using the values in the phi factor file for that realizations - WAT?! I setup this up and tried on a project but it didnt work like I expected (its a lot of fresh code in there and its undocumented so user beware if you try this! haha)

[cnicol-gwlogic]

Cool, thanks. And presumably if I were to not apply a weights ensemble - so I simply use a single set of weights, and phi contributions are defined as the mean prior contributions scaled by the ies_phi_factor_file factors - multimodal should (theoretically) also help with attaining a generally better fit to the obs across the ensemble, because upgrade directions per realz can vary? Whereas without multimodal, directions are defined across the entire ensemble, so one direction - have I got this anywhere near right?

Presumably even without a weights ensemble, the obs-modelled fit should be generally better than if I had not used multimodal, because there is greater flexibility in upgrades for each realisation. Would that be fair to say?

[jtwhite79]

thats a good point - the inverse problem is nonlinear, the multimodal solution should help even without a weight ensemble.

I was running a few experiments with the weight ensemble after this discussion....here is the standard MM solution for that circle problem (just one observation):

So you can see that with nonlinear inverse problem, the MM solution works well even without messing with weights at the realization level. Now lets add another obs, which is just that we also want the final value of parameter 1 to equal its upper bound. Using a standard phi factor file with 50% for each of the two obs we get:

So you can see that the MM solution does exactly what we asked: split the difference between the circle obs and the upper bound of par1. But that's kinda meh, right? So check this out with using a phi factor csv file where 1/3 of the reals are focused on the circle, 1/3 are focused on the upper bound of par1 and 1/3 split the two (and this is without needing to supply a weight ensemble):

That is probably what we are after in practice Im guessing: a posterior that respects more facets of the data...

The settings to make these runs are in this script if you wanna see:

https://github.com/pestpp/pestpp-ies_benchmarks/blob/ea83e259b2f81dcdb5e2d794bba8e33c7703259b/ies_test_part4.py#L1108

and phi factors by reals part is specifically here:

https://github.com/pestpp/pestpp-ies_benchmarks/blob/ea83e259b2f81dcdb5e2d794bba8e33c7703259b/ies_test_part4.py#L1244

Obvs its more complicated in practice, but Im keen to harden this part of ies bc this is usually more what Im looking...so if you wanna try it out, I'll do my best to make it work for you!

[cnicol-gwlogic]

Nice examples! Thanks.

That's a nice point you make re: splitting the difference too. The 50/50 case seems no good in practical terms. I am wondering if in effect this may be similar to splitting the weights across lots of obs groups reasonably evenly. Example - I've got 25-odd obs groups, ranging from 2-5% share each. Behind the scenes, I've taken into consideration lumping these fine groups together into broader groups (eg flows (from duration curves, absolute flows, annual flows, etc) vs heads (at different locations and times of interest) vs drawdowns (similar to heads lumping)). Maybe that is similar to your "split the difference" obs group effect - because I've kind of already split the difference in my weighting strategy - i.e. there are lots of things for pest the upgrades to focus on, and all are weighted similarly.

But I take your point, the 50/50 example scares me a bit - here's hoping I don't end up there!

[jtwhite79]

That example is a bit dramatic bc the two obs are competing...If your model is perfect (!) then in theory you should be able to match all those 25 groups simultaneously...I'm keen to hear how it goes!

[cnicol-gwlogic]

Haha course it's perfect. I'll keep you posted.

[wkitlasten]

Keen to hear more! I too am plagued by "some more outlandish phi components for a few realz" that I suspect are inhibiting my model from attaining a fit commencerate with its level of perfection. (x-axis is log(Tr), so only off by several orders of magnitude, meaning a few thousand years!)

[cnicol-gwlogic]

Hey wes, I guess I was just thinking that the "crazy" / experimental phi factors by realz option could be useful to make the weighting specific to each realz. So for the outlier realz across the ensemble, they would be balanced differently to those that don't have the outliers, and that maybe (maybe) that would be helpful to the upgrade process. But I'm flying blind a little here.

I think you might be able to achieve similar with a weights ensemble anyway. If I were to use the phi factors by realz option, then I would also have to provide a weights ensemble, but I'd make them all the same (eg poor vs good quality data points), and use the phi factors to scale the various obs groups differently across the realz.

Either way, I guess it would mean running the prior first, assessing for outliers, designing that weight and/or phi factors ensemble, and then ploughing on into the upgrade iterations.

[jtwhite79]

I agree with @cnicol-gwlogic - I think a weight ensemble might useful too, Wes. If you could run the prior, then offline setup a weight ensemble that balances each realization so that those outlier realz get the weights they deserve (haha!), then restart with all the req ensembles (and using the multimodal solve), that should help I think. Beware! The multimodal solve is not meant to be used with localization tho!

[wkitlasten]

setup a weight ensemble that balances each realization so that those outlier realz get the weights they deserve (haha!)

The idea being each real will have approx the same phi, but the outliers in the relevant reals will be weighted really low?

The multimodal solve is not meant to be used with localization tho!

Yep, been there!!!

…

-- "Perfect spheres are pointless."

[jtwhite79]

yeah thats the idea - since each realization is being weight balanced individually, those reals with "less than desirable" residuals for some obs will naturally have weights that accommodate those special residuals...at least in theory thats the idea. have you tried a weight ensemble yet?

[wkitlasten]

I have not tried it yet. I need to let it rattle around in my brain a bit more.

…

On Thu, Jun 15, 2023 at 6:23 PM J Dub ***@***.***> wrote: yeah thats the idea - since each realization is being weight balanced individually, those reals with "less than desirable" residuals for some obs will naturally have weights that accommodate those special residuals...at least in theory thats the idea. have you tried a weight ensemble yet? — Reply to this email directly, view it on GitHub <#255 (reply in thread)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ADSJXRGQDFPRGOGMG6UQNC3XLKS5NANCNFSM6AAAAAAZD6P3TA> . You are receiving this because you commented.Message ID: ***@***.***>

-- "Perfect spheres are pointless."

[cnicol-gwlogic]

Just looking at your histograms again, you've got a lot of realz that aren't moving left, not just a few. So maybe what we are talking about won't be of great help there, and I wonder if the issue is:

1. Telling you that your tritium concs can only be explained by those few parameter sets that can get down there. Maybe that's just the way it is. You could try resetting your prior somewhere around the mean of those tritium-sensitive parameters that allow phi to get down to where you want it (?).
2. And/or model structural limitations are hindering those other param sets getting down there (eg. layering, grid resolution, parameter resolution).

Another possibility is that other obs groups are clashing with these tritium ones - and maybe that's where different weighting strategies via multimodal solves might help.

Anyway, sounds an interesting model.