Add files via upload#25
Conversation
Volumetric thermal expansion data
…l expnsion coeffecient NO_JIRA
…s to match format of rest of the repo NO_JIRA
| print("Or, modify this script to match your column names.") | ||
| exit(1) | ||
| # Helper function to process a subset | ||
| def process_and_save(sub_df, output_csv): |
There was a problem hiding this comment.
🚫 [flake8] <302> reported by reviewdog 🐶
expected 2 blank lines, found 0
|
|
||
| # For each T_point, calculate confidence interval | ||
| # Design matrix for quadratic | ||
| X = np.vstack([T_points**2, T_points, np.ones_like(T_points)]).T |
There was a problem hiding this comment.
[flake8] <841> reported by reviewdog 🐶
local variable 'X' is assigned to but never used
| output_df = pd.DataFrame({ | ||
| 'Identifier': np.arange(1, len(T_points) + 1), | ||
| 'Property': [f"Volumetric thermal expansion coefficient @{int(T)}K" for T in T_points], | ||
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, |
There was a problem hiding this comment.
🚫 [flake8] <122> reported by reviewdog 🐶
continuation line missing indentation or outdented
| 'Property': [f"Volumetric thermal expansion coefficient @{int(T)}K" for T in T_points], | ||
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, | ||
| 'Value (1/K)': norm_exp_coeffs, | ||
| # 'Regression Fit (Volume)': volumes_at_points, |
There was a problem hiding this comment.
🚫 [flake8] <122> reported by reviewdog 🐶
continuation line missing indentation or outdented
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, | ||
| 'Value (1/K)': norm_exp_coeffs, | ||
| # 'Regression Fit (Volume)': volumes_at_points, | ||
| # 'Confidence Band (+/-)': conf_band, |
There was a problem hiding this comment.
🚫 [flake8] <122> reported by reviewdog 🐶
continuation line missing indentation or outdented
|
|
||
| # For each T_point, calculate confidence interval | ||
| # Design matrix for quadratic | ||
| X = np.vstack([T_points**2, T_points, np.ones_like(T_points)]).T |
There was a problem hiding this comment.
[flake8] <841> reported by reviewdog 🐶
local variable 'X' is assigned to but never used
| output_df = pd.DataFrame({ | ||
| 'Identifier': np.arange(1, len(T_points) + 1), | ||
| 'Property': [f"Volumetric thermal expansion coefficient @{int(T)}K" for T in T_points], | ||
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, |
There was a problem hiding this comment.
🚫 [flake8] <122> reported by reviewdog 🐶
continuation line missing indentation or outdented
| 'Property': [f"Volumetric thermal expansion coefficient @{int(T)}K" for T in T_points], | ||
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, | ||
| 'Value (1/K)': norm_exp_coeffs, | ||
| # 'Regression Fit (Volume)': volumes_at_points, |
There was a problem hiding this comment.
🚫 [flake8] <122> reported by reviewdog 🐶
continuation line missing indentation or outdented
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, | ||
| 'Value (1/K)': norm_exp_coeffs, | ||
| # 'Regression Fit (Volume)': volumes_at_points, | ||
| # 'Confidence Band (+/-)': conf_band, |
There was a problem hiding this comment.
🚫 [flake8] <122> reported by reviewdog 🐶
continuation line missing indentation or outdented
| output_df.to_csv(output_csv, index=False) | ||
|
|
||
| # Separate SC-XRD and PXRD data | ||
| df[NOTES_COL] = df[NOTES_COL].astype(str) |
There was a problem hiding this comment.
🚫 [flake8] <305> reported by reviewdog 🐶
expected 2 blank lines after class or function definition, found 1
| @@ -0,0 +1,120 @@ | |||
| # import required libraries | |||
| from scipy.stats import t | ||
|
|
||
| # User parameters | ||
| INPUT_CSV = sys.argv[1] # Change as needed |
There was a problem hiding this comment.
Are their going to be instructions on how to run this on the command line? I guess the user won't change this line, rather they are free to specify different filenames which are picked up by the script.
|
|
||
| # User parameters | ||
| INPUT_CSV = sys.argv[1] # Change as needed | ||
| TEMP_COL = 'Temperature/ K' |
There was a problem hiding this comment.
I'm generally not a fan of units in column headers. These can be defined elsewhere, even in a separate 'code book'.
| VOLUME_COL = 'volume / Angstrom-3' | ||
| NOTES_COL = 'notes' | ||
| INCREMENT = 10 # Kelvin | ||
| if "Form" in INPUT_CSV: |
There was a problem hiding this comment.
Form could be a string constant FORM = 'Form', and re-used.
| else: | ||
| name = "" | ||
| # Output file names | ||
| OUTPUT_SCXRD = f'thermal_expansion_coefficient_{name}SCXRD.csv' |
There was a problem hiding this comment.
thermal_expansion_coefficient could be a string constant
| # 'Thermal Expansion Coefficient (dV/dT)': exp_coeffs, | ||
| 'Value (1/K)': norm_exp_coeffs, | ||
| # 'Regression Fit (Volume)': volumes_at_points, | ||
| # 'Confidence Band (+/-)': conf_band, |
There was a problem hiding this comment.
I guess this means the user has to + and - this value against the appropriate predicted value.
| # For each T_point, calculate confidence interval | ||
| # Design matrix for quadratic | ||
| X = np.vstack([T_points**2, T_points, np.ones_like(T_points)]).T | ||
| X_data = np.vstack([temps**2, temps, np.ones_like(temps)]).T |
There was a problem hiding this comment.
On this section of code, quite a few comments.
[1] I think it's trying to compute a 95% confidence interval for a predicted mean because in line 82 there is no +1 inside the sqrt. It looks kind of OK but it constructs things a bit back to front for me (temps**2 as the first column etc in X_data. It'll be fine if everything is consistently ordered).
[2] conf_band looks like a list of numbers rather than intervals. I'd expect the interval to be [y-conf_band, y+conf_band] but I guess the output dataframe only claims a (+/-) here so it's up to the user to do the arithmetic. The script could be changed to do this for them as in actually output [y-conf_band, y+conf_band] .
[3] Unclear what the norm_conf_band is doing but maybe output 'y' is divided by volumes, so the conf interval on y would need to be adjusted to generate a conf. interval on y/volume.
[4] As said, it looks sort of OK but it's doing it all 'by hand' and it would be good to verify against data for which we know the answer. On this point, there is a python model statsmodels. We could try to use this? It seems fairly straighforward to get both confidence and prediction intervals, and fit diagnostics.
Here is some python do fit a quadratic using statsmodels and get confidence and prediction intervals. I am confident about this as I managed to reproduce the exact same confidence and prediction intervals in R.
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.formula.api as smf
# Generate synthetic data
np.random.seed(0)
n = 10
x = np.linspace(0, 10, n)
e = np.random.normal(size=n)
y = 1 + 0.5 * x + 0.8*x**2 + 2 * e
print(f'y {y}')
X = pd.DataFrame({'x':x})
print(f'X {X}')
model = smf.ols(formula = 'y ~ x + np.power(x, 2)', data = {'y':y, 'x':X}).fit()
print(model.summary())
pred = model.get_prediction(X)
pred_summary = pred.summary_frame(alpha=0.05) # 95% confidence intervals
# Extract confidence intervals
ci_lower = pred_summary['mean_ci_lower']
ci_upper = pred_summary['mean_ci_upper']
pi_lower = pred_summary['obs_ci_lower']
pi_upper = pred_summary['obs_ci_upper']
plt.scatter(x, y, label='Data')
plt.plot(x, model.fittedvalues, color='red', label='Fitted Line')
plt.fill_between(x, ci_lower, ci_upper, color='red', alpha=0.3, label='95% CI')
plt.fill_between(x, pi_lower, pi_upper, color='blue', alpha=0.2, label='95% PI')
plt.xlabel('X')
plt.ylabel('Y')
plt.legend()
plt.show()
| if has_scxrd: | ||
| process_and_save(scxrd_df, OUTPUT_SCXRD) | ||
| else: | ||
| print("No SC-XRD data found. Explicitly labelled 'SC-XRD' in the Notes column") |
There was a problem hiding this comment.
Does this print statement text express what you want?
| if has_pxrd: | ||
| process_and_save(pxrd_df, OUTPUT_PXRD) | ||
| else: | ||
| print("No PXRD data found. Explicitly labelled 'PXRD' in the Notes column") |
There was a problem hiding this comment.
Does this print statement text express what you want?
| else: | ||
| print("No PXRD data found. Explicitly labelled 'PXRD' in the Notes column") | ||
| else: | ||
| print("No Explicitly labelled SC-XRD or PXRD data found. Processing all data together.") |
Volumetric thermal expansion data. Making a full PR since this seems to me to be data we can easily add. I assume you want to delete some files in other files though so do that before merging