Comments (16)
wedeling
commented on June 26, 2024
1
Hi @mlchodkowski, thanks for this bug report. I was aware that this can create issues, and some time ago added a quick hack in the SC sampler such that it would actually return integer values. I'll have a look soon
from easyvvuq.
wedeling
commented on June 26, 2024
@mlchodkowski I'm working in this now, sorry it took so long, was on parental leave due to the birth of my son amongst other things. As a first comment, I think the issue is that you can only have one type in a numpy array. So even if chaospy will give a integer distribution, numpy turns it into a float once combined with the other parameters. A quick workaround that should work for the SC sampler is
campaign = uq.Campaign(name="steam-water", params=params, actions=actions,
work_dir=WORK_DIR, verify_all_runs=False)
The verify_all_runs flag just turns off the type check, and since the SC sampler gives integer values this should work. I will try to make a less hacky solution that converts the floats once they are passed from the sample array (sampler.xi_d in the SC case) to the solver. The PCE error you flag is probably due to an incorrectly selected quad rule. Will get back to you later today.
from easyvvuq.
wedeling
commented on June 26, 2024
@mlchodkowski I think it should be fixed now, without the workaround mentioned above. If you pull the issue401 branch (see #402) you can check. If so, I'll merge it into the dev branch.
from easyvvuq.
wedeling
commented on June 26, 2024
@mlchodkowski I went ahead and merged the pull request. You can now check your implementation by just pulling the latest dev. If there are still issues let me know.
from easyvvuq.
mlchodkowski
commented on June 26, 2024
Hi @wedeling,
first of all, thanks for handling this problem. I've successfully pulled the changes into my local repo. AFAIK the SC sampler now works correctly, but I'm still testing the changes on my models.
Right now I've encountered a problem with PCE sampler. AFAIK the invalid quadrature rule that I probably use should now be automatically replaced when using DiscreteUniform distribution. However I've encountered a bug like this:
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/easyvvuq/sampling/pce.py", line 257, in __init__
self._nodes, self._weights = cp.generate_quadrature(order=polynomial_order,
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/frontend.py", line 135, in generate_quadrature
return sparse_grid(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/sparse_grid.py", line 81, in sparse_grid
x_lookup, w_lookup = _construct_lookup(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/sparse_grid.py", line 146, in _construct_lookup
(abscissas,), weights = chaospy.generate_quadrature(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/frontend.py", line 152, in generate_quadrature
abscissas, weights = _generate_quadrature(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/frontend.py", line 248, in _generate_quadrature
abscissas, weights = quad_function(order, dist, **parameters)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/discrete.py", line 48, in discrete
order = numpy.where(growth, numpy.where(order > 0, 2**order, 0), order)
ValueError: invalid literal for int() with base 10: 'false'
Additional context (this is my internal JSON format, which is further casted to a proper vary dict):
"sampled_params": [
{
"name": "initial_dose",
"sampling_space": {
"param_type": "integer",
"default": 1000
},
"distribution": {
"dist_type": "discrete_uniform",
"lower": 100,
"upper": 1000
}
},
{
"name": "k12",
"sampling_space": {
"param_type": "float",
"default": 0.2
},
"distribution": {
"dist_type": "uniform",
"lower": 0.3,
"upper": 0.7
}
},
{
"name": "k21",
"sampling_space": {
"param_type": "float",
"default": 0.1
},
"distribution": {
"dist_type": "uniform",
"lower": 0.3,
"upper": 0.7
}
},
{
"name": "k10",
"sampling_space": {
"param_type": "float",
"default": 0.15
},
"distribution": {
"dist_type": "uniform",
"lower": 0.3,
"upper": 0.7
}
},
{
"name": "time_steps",
"sampling_space": {
"param_type": "integer",
"default": 100
}
},
{
"name": "time_interval",
"sampling_space": {
"param_type": "integer",
"default": 1
}
}
],
"method": {
"method_type": "pce",
"polynomial_order": "3"
}
I've verified that the problem is not in this format, nor the casting-to-vary step.
from easyvvuq.
wedeling
commented on June 26, 2024
I've reopened the issue. Are you using a sparse grid or a standard tensor grid?
from easyvvuq.
mlchodkowski
commented on June 26, 2024
I don't use sparse grid (default arg sparse=False for PCE)
from easyvvuq.
wedeling
commented on June 26, 2024
It's a bit hard to see from here, could you post your vary dict and the polynomial order you use?
from easyvvuq.
mlchodkowski
commented on June 26, 2024
Sure, here it is:
params = {
"initial_dose": {"type": "integer", "default": 100},
"k12": {"type": "float", "default": 0.2},
"k21": {"type": "float", "default": 0.1},
"k10": {"type": "float", "default": 0.15},
"time_steps": {"type": "integer", "default": 100},
"time_interval": {"type": "integer", "default": 1},
}
vary = {
"initial_dose": cp.DiscreteUniform(100, 1000),
"k12": cp.Uniform(0.3, 0.7),
"k21": cp.Uniform(0.3, 0.7),
"k10": cp.Uniform(0.3, 0.7)
}
campaign.set_sampler(uq.sampling.PCESampler(vary=vary, polynomial_order=3))
EDIT: And here is the model that I'm using:
#!/usr/bin/env python3
import dataclasses
import json
import random
import sys
import numpy as np
import pandas as pd
@dataclasses.dataclass
class Config:
# Parameters
initial_dose: int # Initial drug dose
k12: float # Rate constant from central to peripheral compartment
k21: float # Rate constant from peripheral to central compartment
k10: float # Elimination rate constant
time_steps: int # Number of time steps
time_interval: int # Time interval between steps
@staticmethod
def decode(cfg):
return Config(**cfg)
def stochastic_term(amp):
return amp * random.uniform(0.1, 1)
if __name__ == "__main__":
if len(sys.argv) < 2:
print("ERROR: Config file not provided")
sys.exit(1)
with open(sys.argv[1]) as fh:
config = json.load(fh, object_hook=Config.decode)
amp = 0.15
# Initialize concentration arrays for central and peripheral compartments
central_concentration = np.zeros(config.time_steps)
peripheral_concentration = np.zeros(config.time_steps)
central_concentration[0] = config.initial_dose
# Simulation loop
for t in range(1, config.time_steps):
# Calculate rate of change for each compartment
dC1_dt = config.k21 * peripheral_concentration[t - 1] - config.k12 * central_concentration[
t - 1
] * stochastic_term(amp)
dC2_dt = (
config.k12 * central_concentration[t - 1]
- config.k21 * peripheral_concentration[t - 1]
- config.k10 * peripheral_concentration[t - 1]
) * stochastic_term(amp)
# Update concentrations using Euler's method
central_concentration[t] = central_concentration[t - 1] + dC1_dt * config.time_interval
peripheral_concentration[t] = peripheral_concentration[t - 1] + dC2_dt * config.time_interval
# Plot the results
time = np.arange(0, config.time_steps * config.time_interval, config.time_interval)
df = pd.DataFrame({"time": time, "central": central_concentration, "peripheral": peripheral_concentration})
df.to_csv("output.csv", index=False)
from easyvvuq.
mlchodkowski
commented on June 26, 2024
@wedeling Hi, I wanted to inform you that I've encountered the same bug also with SC sampler. The vary is the same as above, sampler is SCSampler with polynomial_order=3. I've tested this on you latest commit (0d27a63) to 'dev' branch.
Error:
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/easyvvuq/sampling/stochastic_collocation.py", line 106, in __init__
self.check_max_quad_level()
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/easyvvuq/sampling/stochastic_collocation.py", line 247, in check_max_quad_level
xi_i, wi_i = cp.generate_quadrature(order + j,
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/frontend.py", line 152, in generate_quadrature
abscissas, weights = _generate_quadrature(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/frontend.py", line 248, in _generate_quadrature
abscissas, weights = quad_function(order, dist, **parameters)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/quadrature/clenshaw_curtis.py", line 66, in clenshaw_curtis
order = numpy.where(growth, numpy.where(order > 0, 2**order, 0), order)
ValueError: invalid literal for int() with base 10: 'false'
from easyvvuq.
wedeling
commented on June 26, 2024
OK, will have a look today.
from easyvvuq.
wedeling
commented on June 26, 2024
It seems to work on my end, I get results like this (1st-order Sobol indices):
I used the following script, when storing your model as model.py in the same directory
import os
import easyvvuq as uq
import chaospy as cp
import matplotlib.pyplot as plt
from easyvvuq.actions import CreateRunDirectory, Encode, ExecuteLocal, Decode, Actions
WORK_DIR = '/tmp'
OUT_COLS = ["peripheral","central"]
params = {
"initial_dose": {"type": "integer", "default": 100},
"k12": {"type": "float", "default": 0.2},
"k21": {"type": "float", "default": 0.1},
"k10": {"type": "float", "default": 0.15},
"time_steps": {"type": "integer", "default": 100},
"time_interval": {"type": "integer", "default": 1},
}
vary = {
"initial_dose": cp.DiscreteUniform(100, 1000),
"k12": cp.Uniform(0.3, 0.7),
"k21": cp.Uniform(0.3, 0.7),
"k10": cp.Uniform(0.3, 0.7)
}
encoder = uq.encoders.GenericEncoder(
template_fname="config.template", delimiter="$", target_filename="input.json")
decoder = uq.decoders.SimpleCSV(target_filename="output.csv", output_columns=OUT_COLS)
execute = ExecuteLocal(f"python3 {os.getcwd()}/model.py input.json")
actions = Actions(CreateRunDirectory(WORK_DIR, flatten=True), Encode(encoder), execute, Decode(decoder))
campaign = uq.Campaign(name="test", params=params, actions=actions, work_dir=WORK_DIR)
# sampler = uq.sampling.SCSampler(vary=vary, polynomial_order=3)
sampler = uq.sampling.PCESampler(vary=vary, polynomial_order=3)
campaign.set_sampler(sampler)
campaign.execute().collate()
data_frame = campaign.get_collation_result()
# analysis = uq.analysis.SCAnalysis(sampler=sampler, qoi_cols=OUT_COLS)
analysis = uq.analysis.PCEAnalysis(sampler=sampler, qoi_cols=OUT_COLS)
campaign.apply_analysis(analysis)
results = campaign.get_last_analysis()
sobols = results.sobols_first()
fig = plt.figure(figsize=[8, 4])
ax1 = fig.add_subplot(121, title=OUT_COLS[0], xlabel="time", ylabel="1st order Sobol index")
ax2 = fig.add_subplot(122, title=OUT_COLS[1], xlabel="time", sharey=ax1)
for param in vary.keys():
ax1.plot(sobols[OUT_COLS[0]][param], label=param)
ax2.plot(sobols[OUT_COLS[1]][param], label=param)
plt.legend(loc=0)
plt.tight_layout()
plt.show()
The config.template file is as follows, also stored in the same directory:
{
"initial_dose": $initial_dose,
"k12": $k12,
"k21": $k21,
"k10": $k10,
"time_steps": $time_steps,
"time_interval": $time_interval
}
Could you perhaps run this, and pinpoint at what line it goes wrong?
from easyvvuq.
mlchodkowski
commented on June 26, 2024
Hi, I've recreated my virtualenv in which i installed the newest easyvvuq. I've run your code from above and I got this error in campaign.apply_analysis(analysis):
/mnt/storage_2/project_data/grant_571/muqsa-app/.venv/lib/python3.10/site-packages/chaospy/descriptives/correlation/pearson.py:46: RuntimeWarning: invalid value encountered in divide
return numpy.where(vvar > 0, cov/vvar, 0)
Error LinAlgError for peripheral when computing cp.QoI_Dist()
from easyvvuq.
mlchodkowski
commented on June 26, 2024
@wedeling I have succesfully localized the bug. The original error was caused by wrong type set in my my pydantic's models for my internal json file. I'm really sorry for the trouble. The error is now fixed. The implementation of SC and PCE samplers allowed for passing the boolean parameters as str's ('false' and 'true' instead of False and True). This was causing the invalid literal for int() with base 10: 'false' on my side.
However the above error (Error LinAlgError) still persist in some cases + there is a new one.
/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/chaospy/descriptives/correlation/pearson.py:46: RuntimeWarning: invalid value encountered in divide
return numpy.where(vvar > 0, cov/vvar, 0)
/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/muqsa/runner.py:557: FutureWarning: Setting an item of incompatible dtype is deprecated and will raise in a future error of pandas. Value '' has dtype incompatible with float64, please explicitly cast to a compatible dtype first.
result.fillna("", inplace=True)
Traceback (most recent call last):
File "/opt/exp_soft/local/skylake/python/3.10.7/lib/python3.10/runpy.py", line 196, in _run_module_as_main
return _run_code(code, main_globals, None,
File "/opt/exp_soft/local/skylake/python/3.10.7/lib/python3.10/runpy.py", line 86, in _run_code
exec(code, run_globals)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/muqsa/runner.py", line 580, in <module>
r = Runner(sys.argv[1]).run()
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/muqsa/runner.py", line 64, in run
self.analyse(campaign, decoder.output_columns, vary, um.method)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/muqsa/runner.py", line 365, in analyse
Runner.analyse_in_campaign(campaign, qoi_columns, method.__root__.method_type)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/muqsa/runner.py", line 572, in analyse_in_campaign
results.plot_sobols_treemap(q, filename=f"results/{q}_sobols_treemap_{postfix}.svg")
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/easyvvuq/analysis/results.py", line 464, in plot_sobols_treemap
squarify.plot(sizes=values, label=keys, color=colors, ax=ax, pad=True)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 244, in plot
rects = padded_squarify(normed, 0, 0, norm_x, norm_y)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 147, in padded_squarify
rects = squarify(sizes, x, y, dx, dy)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 136, in squarify
return layout(current, x, y, dx, dy) + squarify(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 128, in squarify
while i < len(sizes) and worst_ratio(sizes[:i], x, y, dx, dy) >= worst_ratio(
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 86, in worst_ratio
for rect in layout(sizes, x, y, dx, dy)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 48, in layout
layoutrow(sizes, x, y, dx, dy) if dx >= dy else layoutcol(sizes, x, y, dx, dy)
File "/home/users/mchodkowski/grant_571/project_data/muqsa-app/.venv/lib64/python3.10/site-packages/squarify/__init__.py", line 38, in layoutcol
height = covered_area / dx
ZeroDivisionError: float division by zero
I suspect the issue is in the results that my model generates some bad or incompatible results. Maybe it's worth to handle this exception?
from easyvvuq.
wedeling
commented on June 26, 2024
@mlchodkowski no worries, glad you got it working. The LinAlgError seems to be related to the Sobol indices. Perhaps the Sobol indices are computed somewhere where there is zero variance (perhaps a boundary condition or something?). In this case you would get a divide by zero. I could look into catching this exception.
from easyvvuq.
wedeling
commented on June 26, 2024
Closing this
from easyvvuq.
Related Issues (20)
- Failures when running the non-regression tests manually (1)
- Failures when running the non-regression tests manually (2)
- Should we be getting 1754 warnings when running the non-regression tests manually?
- Possible issue when running test script (basic_tutorial_qcgpj.py) in tutorials
- Fix test_mcmc HOT 1
- Update testing infrastructure?
- Fix test_surrogate_workflow HOT 2
- Fix triggering of LGTM / travis-ci tests
- help needed in converting a hack to a proper implementation for coupling easyvvuq and muscle3
- PDoc fails with cryptic error message when misformatting Returns field in API documentation
- Code Cleanup HOT 1
- Excessive amounts of memory are sometimes used by the PCE analysis phase HOT 1
- Fix code scanning alert - Wrong name for an argument in a call
- Fix code scanning alert - Wrong name for an argument in a call
- Fix code scanning alert - Jinja2 templating with autoescape=False HOT 2
- Fix code scanning alert - Potentially uninitialized local variable
- Fix code scanning alert - Potentially uninitialized local variable
- assertion failure while running easyvvuq_fusion_tutorial.py HOT 1
- Unhandled divide by zero case in sc_analysis.py
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google ❤️ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.
from easyvvuq.