Tibia shaft (Harden 2021)

Tibia shaft (Harden 2021)#

Validation model information#

Performed by: Elena Weissenbacher, Nico Erlinger
Reviewed by: Corina Klug

Added to VIVA+ Validation Catalog on: 2022-11-31

Last modified: 2023-11-23

Model version (this notebook run for): 0.3.2

The Jupyter notebooks are licensed under Creative Commons Attribution 4.0 International License CCBY

Experiments by Harden et al. (2021)#

Harden, A.; Kang, Y. S.; Hunter, R.; Bendig, A.; Bolte, J.; Eckstein, N.; Smith, A.; Agnew, A. 2021. Preliminary Sex-specific Relationships between Peak Force and Cortical Bone Morphometrics in Human Tibiae Subjected to Lateral Loading. In Proccedings of the IRCOBI Conference 2021, Online. http://www.ircobi.org/wordpress/downloads/irc21/pdf-files/2133.pdf

Summary#

Five female and five male isolated tibiae were subjeted to four-point bending.

Information on the subjects/specimens#

Test number	Sex	Donor age [yr]	Mechanical span [mm]
Tib-001	M	64	191
Tib-002	M	63	188
Tib-003	M	77	166
Tib-004	F	84	161
Tib-005	F	86	161
Tib-006	F	102	160
Tib-007	M	73	165
Tib-008	M	74	183
Tib-009	F	85	162
Tib-010	F	89	148

Loading and Boundary Conditions#

A constant impactor velocity of 6 m/s was prescribed. Mechanical span of the potting was adjusted to experimental data for each simulation.

Experimental responses#

The force-time history plots published by Harden et al. (2021) were digitised with WebPlotDigitizer v4.4 (https://automeris.io/WebPlotDigitizer).

Other notes for simulation#

The simulation setups will be shared as soon as approved by the authors of the study.

processed_data_dir = f'data/processed/2022-11-30_Name_of_person'

name ='Name_of_person'
date = datetime.date.today().strftime('%Y-%m-%d')
dynasaur_output_file_name = 'Dynasaur_output.csv'

##Results

Impactor force plots#

../_images/ffeae8eb13a62ff0f3ffd950727818d541bdf66368ec52c5a254505cca74fa01.png

Plot strain-time histories#

../_images/aefa4c973d6b2ed2a34dcb38160342eba11b6e97969d6875e27ddad44f0668f6.png

Engergy-time plots#

../_images/51e5ad667e08aa4cde52f439f1cb7f07ecfb4456b80d0a65aaa554a9054c7d2a.png

ISO18571 objective rating for force-time histories#

Show code cell source Hide code cell source

simulation_list = ['Tib-001','Tib-002','Tib-003','Tib-004','Tib-005','Tib-006','Tib-007','Tib-008','Tib-009','Tib-010']

def calculate_iso_score(processed_data_dir, experiment_data_file, sim_x_data, sim_y_data, sim_data_for_calculation_y, exp_x_data, exp_y_data):
    ratings = []
    for i in simulation_list:
        experimental_data = pd.read_csv(experiment_data_file, delimiter=';', header=[0,1,2,3], decimal='.')
        processed_data_path = os.path.join(processed_data_dir, i).replace('\\', '/')
        simData = pd.read_csv(os.path.join(processed_data_path, dynasaur_output_file_name), 
                                    delimiter=';', na_values='-', header = [0,1,2,3]) 
              
        time_ref = np.array(experimental_data[(i,) + exp_x_data]).flatten()
        x_ref = np.array(experimental_data[(i,) + exp_y_data]).flatten()

        time_comp=np.array(simData[sim_x_data]).flatten()
        x_comp = simData[sim_y_data]
        
        if sim_data_for_calculation_y is not None:
            x_comp = x_comp + simData[sim_data_for_calculation_y]
        x_comp=np.array(x_comp).flatten() 
        
        time_comp = np.delete(time_comp,np.s_[150:])
        x_comp = np.delete(x_comp, np.s_[150:])

        ind_nan = []
        ind_nan = np.where(np.isnan(x_ref))
        
        if len(ind_nan[0]) > 0:
            x_ref = x_ref[0:ind_nan[0][0]-1]
            x_time_ref = time_ref[0:ind_nan[0][0]-1]
            x_comp = x_comp[0:ind_nan[0][0]-1]
            x_time_comp = time_comp[0:ind_nan[0][0]-1]
        else:
            x_time_ref = time_ref
            x_time_comp = time_comp
            
        if len(x_ref) > len(x_comp):
            ind_end = len(x_comp)
            x_ref = x_ref[0:ind_end]
            x_time_ref = x_time_ref[0:ind_end]
            
        elif len(x_comp) > len(x_ref):
            ind_end = len(x_ref)
            x_comp = x_comp[0:ind_end]
            x_time_comp = x_time_comp[0:ind_end]   

        ref = np.vstack((x_time_ref, x_ref)).T
        comp = np.vstack((x_time_comp, x_comp)).T

        iso_rating = ISO18571(reference_curve=ref, comparison_curve=comp)
        
        ratings.append([i,
        iso_rating.corridor_rating(), 
        iso_rating.phase_rating(),
        iso_rating.magnitude_rating(),
        iso_rating.slope_rating(),
        iso_rating.overall_rating()])

    rating_array = []
    rating_array = np.append(rating_array, ratings)
    rating_array = rating_array.transpose()
    rating_array = np.reshape(rating_array, (10,6))

    return(rating_array)

iso_ratings = calculate_iso_score(
sim_x_data = ('TIBIA', 'Contact_Force_proximal_z', 'time'),
sim_y_data = ('TIBIA', 'Contact_Force_proximal_z', 'force'),
sim_data_for_calculation_y = ('TIBIA', 'Contact_Force_distal_z', 'force'),
experiment_data_file = 'data/experiment/Harden_2021_testdata.csv',
exp_x_data = ('time', 'ms'),
exp_y_data = ('force', 'kN'),
processed_data_dir = f'data/processed/2022-11-30_Name_of_person'
)

df = pd.DataFrame(iso_ratings, columns=['Simulation', 'Corridor', 'Phase', 'Magnitude', 'Slope', 'Overall'])
display(df)

	Simulation	Corridor	Phase	Magnitude	Slope	Overall
0	Tib-001	0.678	0.9	0.742	0.34	0.668
1	Tib-002	0.637	0.767	0.635	0.676	0.67
2	Tib-003	0.714	0.833	0.898	0.599	0.752
3	Tib-004	0.722	0.9	0.718	0.598	0.732
4	Tib-005	0.757	0.833	0.824	0.607	0.756
5	Tib-006	0.69	0.8	0.813	0.532	0.705
6	Tib-007	0.646	0.9	0.741	0.655	0.718
7	Tib-008	0.753	0.967	0.91	0.705	0.818
8	Tib-009	0.601	0.733	0.712	0.491	0.628
9	Tib-010	0.63	0.833	0.432	0.573	0.619

Quasistatic simulations#