rmaia / pavo

tools for the analysis of color data in R
http://pavo.colrverse.com
GNU General Public License v2.0
68 stars 17 forks source link

Issue with how photopigment absorbance is processed in sensmodel #256

Open hneyster opened 2 months ago

hneyster commented 2 months ago

Thanks so much for all your work on this amazing package that has been so useful to so many!

Colleagues and I have identified what seems to be an issue with how the photopigment absorbance is processed in `sensmodel.R``:

The Govardovskii et al. (2000) equations are used to define the # Sensitivities w/o oil droplets and denoted peak in lines 94-100 in sensmodel.R. These equations provide estimates of photopigement ABSORBANCE [i.e., the log ratio of incident to transmitted light] (see Govardovskii et al., 2000 and Endler & Mielke, 2005).

In order to integrate photo pigment ABSORBANCE with oil droplet transmittance, one has to first convert the photo pigment ABSORBANCE to ABSORPTANCE [i.e., 1-transmittance, on fractional rather than log scale]. Otherwise one is multiplying a log-ratio by a fractional ratio.

Endler & Mielke (2005) show how to carry out this transformation in their equation 3, where $1-10^{-ad_rG_r(\lambda)}$ describes how the Govardovskii log-ratio ABSORBANCE ( $G_r(\lambda)$ ) is converted to fractional-ratio ABSORPTANCE.

The resulting ABSORPTANCE can then be multiplied by the transmittance of the oil droplet etc to obtain the photoreceptor capture probability.

However, in sensmodel.R line 128, the Govardovskii ABSORBANCE is multiplied directly with the oil droplet transmittance.

peak <- peak * T.oil

This leads to an incorrect estimate of photoreceptor sensitivity.

Here's some code to calculate and plot how ABSORBANCE and ABSORPTANCE have slightly different shapes:

# setting some parameters: 
peaksens = 570 # lambda_max of bird photoreceptor 1 
peaksens2 =412 # lambda_max of bird photorecepto 2
wl <- 300:700 # wavelengths

# eye parameters: increasing the size of these parameters increases the difference between ABSORBANCE and ABSORPTANCE
alpha <- .014 # pigment specific absorbance; see Endler & Mielke 2005; representative value from Hart & Vorobyev 2005
d_r <- 16 # outer segment path length from Endler & Mielke 2005; representative value from Hart & Vorobyev 2005 

## Absorpbance: 
peak1_ab <- 1/(exp(69.7 * (0.8795 + 0.0459 * exp(-(peaksens - 
           wl[1])^2/11940) - (peaksens/wl))) + exp(28 * 
          (0.922 - peaksens/wl)) + exp(-14.9 * (1.104 - 
            (peaksens/wl))) + 0.674) +
  0.26 * exp(-((wl - (189 + 0.315 * peaksens))/(-40.5 + 
     0.195 * peaksens))^2) # govardovskii pigment absorbance 

peak2_ab <- 1/(exp(69.7 * (0.8795 + 0.0459 * exp(-(peaksens2 - 
 wl[1])^2/11940) - (peaksens2/wl))) + exp(28 * 
(0.922 - peaksens2/wl)) + exp(-14.9 * (1.104 - 
  (peaksens2/wl))) + 0.674) + 
  0.26 * exp(-((wl - (189 + 0.315 * peaksens2))/(-40.5 + 
       0.195 * peaksens2))^2)

# Absorptance: 

peak1_ap <- (1-10^(-alpha*d_r*peak1_ab)) 
peak2_ap <- (1-10^(-alpha*d_r*peak2_ab))

op <- par(mar=c(5, 6, 4, 2) + 0.1)
plot(wl,peak2_ap/max(peak2_ap),ylim=c(0,1.1),ylab = 'normalized absorptance (black) \n normalized absorbance(red)')
points(wl,peak1_ap/max(peak1_ap))
points(wl,peak2_ab/max(peak2_ab), col='red')
points(wl,peak1_ab/max(peak1_ab),col='red')

Created on 2024-06-17 with reprex v2.1.0

It also looks like there are a couple of places in the documentation and axes labels that are erroneously labeled 'absorbance.' Two examples: Figure 3.4 in the documentation is labeled with 'absorbance' on the y-axis, but is usually called capture probability or spectral sensitivity etc. (see, e.g., Endler & Mielke 2005; Hart & Vorobyev 2005).

Section 3.1.4.4 also erroneously uses 'absorbance' in this sentence:

...columns being the absorbance at each wavelength for each cone type (see below for an example).

Absorbance and absorptance can be tricky to get right, but to both be mathematically and biologically correct, I recommend converting photo pigment absorbance to absorptance, as described by Endler & Mielke (2005), before multiplying by oil droplet transmittance in sensmodel and double-checking the use of 'absorbance' in the documentation and axes labels.

References:

Endler John A. And Mielke, P. W. (2005). Comparing entire colour patterns as birds see them. Biological Journal of the Linnean Society, 86(4), 405–431. https://doi.org/10.1111/j.1095-8312.2005.00540.x

Govardovskii, V. I., Fyhrquist, N., Reuter, T., Kuzmin, D. G., & Donner, K. (2000). In search of the visual pigment template. Visual Neuroscience, 17(4), 509–528. https://doi.org/10.1017/s0952523800174036

Hart, N. S., & Vorobyev, M. (2005). Modelling oil droplet absorption spectra and spectral sensitivities of bird cone photoreceptors. Journal of Comparative Physiology A, 191(4), 381–392. https://doi.org/10.1007/s00359-004-0595-3

Bisaloo commented 1 month ago

Thanks a lot for the detailed report!

I went through Endler & Mielke (2005) methods again and your report makes sense to me. I'll try to tackle it in the coming weeks in a larger refactor of sensmodel().

Feel free to also send a PR if you would like a faster fix for the issue at hand, without having to wait for the refactor :relaxed:. Thanks!