for func in el_reaction.rate_laws[0].expression.functions:                
                for species in func.expression.species:
                    el_species_counts[species.id] = species.distribution_init_concentration.mean
                for obs in func.expression.observables:
                    for species in obs.expression.species: