Shape transition error
40 views
Skip to first unread message
Shane Morris
Apr 14, 2026, 9:37:54 AMApr 14
to NicheMapR
Hi Mike,
I've noticed an error that occurs when an endotherm (as an ellipsoid) transitions from a shape (SHAPE_B) where their percentage of subcutaneous fat is enough to have a non-zero fat thickness (output FAT_THICK) to a shape where their fat thickness becomes 0. It does not occur when subcutaneous fat is not present (SUBQFAT=0). It appears to be an issue with the calculation of the flesh conductivity (K_FLESH). The flesh conductivity does not change when the SHAPE_B parameter changes. Could be something to do with how the calculation of whether FAT_THICK is zero or non-zero occurs in the GEOM_ENDO function which is after the skin conductivity value (AK1_LAST) is calculated in SOLVENDO.f? Below is the code to explain what I mean. It is for a hibernating brown bear.
Cheers,
Shane.
# load libraries
library(NicheMapR); library(tidyverse)
# Input variables for a hibernating bear
TAs <- seq(-30, 30, 1)
AMASS <- 135
Q10_torpor <- 1
QBASAL <- 0.3509*AMASS
TC <- 37.5
TC_MAX <- 40
TC_MIN <- 30
TC_INC <- 0.01
# The variables the issue is related to
FATPCT <- 20
SUBQFAT <- 1
SHAPE_B <- 1.3
SHAPE_B_MAX <- 6
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out1 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out1$TAs <- TAs
# Plot the metabolic rate
ggplot(endo.out1)+ geom_line(aes(x=TAs, y=enbal.QGEN))+ theme_classic()
# Why this peak at 7 degrees?
ggplot(endo.out1)+ geom_line(aes(x=TAs, y=treg.SHAPE_B))+ theme_classic()
ggplot(endo.out1)+ geom_line(aes(x=TAs, y=morph.AREA))+ theme_classic()
# Also why does the surface area decrease when posture changes? The surface area of a more spherical object should be lower than one that is more ellipsoidal. Is it to do with the fur flattening, thus decreasing the surface area?
# This blip in metabolism is tied to subcutaneous fat. Here I will plot a model with no subcutaneous fat, one with 50% extra fat and one with 50% less fat
SUBQFAT <- 0
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out2 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out2$TAs <- TAs
# This blip is clearly tied to the fat percentage. It dramatically rises with more fat and dramatically decreases with
SUBQFAT <- 1
FATPCT <- 30
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out3 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out3$TAs <- TAs
FATPCT <- 10
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out4 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out4$TAs <- TAs
endo.out2 %>% mutate(FAT="None")
# Plot the metabolic rate
(endo.out1 %>% mutate(FAT="Normal") %>% bind_rows(endo.out2 %>% mutate(FAT="None"), endo.out3 %>% mutate(FAT="50% extra"), endo.out4 %>% mutate(FAT="50% less"))) %>%
ggplot()+ geom_line(aes(x=TAs, y=enbal.QGEN, group=FAT, colour=FAT), linewidth=2)+ theme_classic()
# It appears to be due to the transition from ellipsoid ratios 4.5 to 4.6
# If I go back to the original fat variables and keep the shape always at 4 or 5 it does not happen
FATPCT <- 20.1
SHAPE_B <- 4.5
SHAPE_B_MAX <- 4.5
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out5 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out5$TAs <- TAs
SHAPE_B <- 4.5
SHAPE_B_MAX <- 4.6
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out6 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out6$TAs <- TAs
SHAPE_B <- 4.6
SHAPE_B_MAX <- 4.6
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out7 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out7$TAs <- TAs
(endo.out5 %>% mutate(SHAPE="Just 4.5") %>% bind_rows(endo.out6 %>% mutate(SHAPE="4.5 to 4.6"), endo.out7 %>% mutate(SHAPE="Just 4.6"))) %>%
ggplot()+ geom_line(aes(x=TAs, y=enbal.QGEN, group=SHAPE, colour=SHAPE), linewidth=2)+ theme_classic()
# Plotting all variables to look at patterns that may explain
for(i in 1:ncol(endo.out5)){
min5 <- min(endo.out5[, i], na.rm=TRUE)
max5 <- max(endo.out5[, i], na.rm=TRUE)
plot(TAs, endo.out5[, i], type="l", main=paste0("Variable ", colnames(endo.out5)[i]), ylim=c(min5*0.5, max5*1.5))
lines(TAs, endo.out6[, i], col="blue")
lines(TAs, endo.out7[, i], col="red")
abline(v=-4, col="grey")
abline(v=11, col="grey")
}
# From looking at the variables. The only one I can see that is a candidate is the K_FLESH
(endo.out5 %>% mutate(SHAPE="Just 4.5") %>% bind_rows(endo.out6 %>% mutate(SHAPE="4.5 to 4.6"), endo.out7 %>% mutate(SHAPE="Just 4.6"))) %>%
ggplot()+ geom_line(aes(x=TAs, y=treg.K_FLESH, group=SHAPE, colour=SHAPE), linewidth=2)+ theme_classic() +
geom_vline(xintercept=-4, colour="grey")+ geom_vline(xintercept=11, colour="grey")
# it appears to be caused by the delay in the increase in flesh conductivity
# Can see here that the transition in SHAPE_B occurs but the skin conductivity doesn't
endo.out6 %>%
ggplot()+ geom_line(aes(x=TAs, y=treg.SHAPE_B), linewidth=2)+ theme_classic() +
geom_vline(xintercept=-4, colour="grey")+ geom_vline(xintercept=11, colour="grey")
library(NicheMapR); library(tidyverse)
# Input variables for a hibernating bear
TAs <- seq(-30, 30, 1)
AMASS <- 135
Q10_torpor <- 1
QBASAL <- 0.3509*AMASS
TC <- 37.5
TC_MAX <- 40
TC_MIN <- 30
TC_INC <- 0.01
# The variables the issue is related to
FATPCT <- 20
SUBQFAT <- 1
SHAPE_B <- 1.3
SHAPE_B_MAX <- 6
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out1 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out1$TAs <- TAs
# Plot the metabolic rate
ggplot(endo.out1)+ geom_line(aes(x=TAs, y=enbal.QGEN))+ theme_classic()
# Why this peak at 7 degrees?
ggplot(endo.out1)+ geom_line(aes(x=TAs, y=treg.SHAPE_B))+ theme_classic()
ggplot(endo.out1)+ geom_line(aes(x=TAs, y=morph.AREA))+ theme_classic()
# Also why does the surface area decrease when posture changes? The surface area of a more spherical object should be lower than one that is more ellipsoidal. Is it to do with the fur flattening, thus decreasing the surface area?
# This blip in metabolism is tied to subcutaneous fat. Here I will plot a model with no subcutaneous fat, one with 50% extra fat and one with 50% less fat
SUBQFAT <- 0
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out2 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out2$TAs <- TAs
# This blip is clearly tied to the fat percentage. It dramatically rises with more fat and dramatically decreases with
SUBQFAT <- 1
FATPCT <- 30
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out3 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out3$TAs <- TAs
FATPCT <- 10
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out4 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out4$TAs <- TAs
endo.out2 %>% mutate(FAT="None")
# Plot the metabolic rate
(endo.out1 %>% mutate(FAT="Normal") %>% bind_rows(endo.out2 %>% mutate(FAT="None"), endo.out3 %>% mutate(FAT="50% extra"), endo.out4 %>% mutate(FAT="50% less"))) %>%
ggplot()+ geom_line(aes(x=TAs, y=enbal.QGEN, group=FAT, colour=FAT), linewidth=2)+ theme_classic()
# It appears to be due to the transition from ellipsoid ratios 4.5 to 4.6
# If I go back to the original fat variables and keep the shape always at 4 or 5 it does not happen
FATPCT <- 20.1
SHAPE_B <- 4.5
SHAPE_B_MAX <- 4.5
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out5 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out5$TAs <- TAs
SHAPE_B <- 4.5
SHAPE_B_MAX <- 4.6
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out6 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out6$TAs <- TAs
SHAPE_B <- 4.6
SHAPE_B_MAX <- 4.6
ptm <- proc.time()
endo.out <- lapply(1:length(TAs), function(x) {
endoR(TA = TAs[x], VEL = 0.1, RH = 5, QSOLR = 0,
AMASS = AMASS, Q10 = Q10_torpor, QBASAL = QBASAL,
TC = TC, TC_MAX = TC_MAX, TC_MIN = TC_MIN,
SHAPE_B = SHAPE_B,SHAPE_B_MAX = SHAPE_B_MAX,
FATPCT = FATPCT, SUBQFAT = SUBQFAT,
TORPOR = 1,
SHAPE=4, SAMODE=2, DHAIRD = 0.000076, DHAIRV = 0.000076, LHAIRD = 0.0759,
LHAIRV = 0.0828,ZFURD = 0.0559,ZFURV = 0.0628, RHOD = 8.5e+6, RHOV = 8.5e+6*0.75,
REFLD = 0.430, REFLV = 0.29, PZFUR = 0.1, PVEN = 0.25,UNCURL = 0.1,
PCTWET = 0.1, PCTWET_MAX = 10,PCTWET_INC = 0.25, DELTAR = 5,EXTREF = 25,
PANT_INC = 0.1, PANT_MAX = 10,PANT_MULT = 1,AK1 = 1.5,AK2 = 0.23,AK1_MAX = 2.8,
AK1_INC = 0.1,ANDENS = 1033.6, RQ = 0.7 )
})
proc.time() - ptm
# combine data frame and add ambient temperature
endo.out7 <- do.call("rbind", lapply(endo.out, data.frame))
endo.out7$TAs <- TAs
(endo.out5 %>% mutate(SHAPE="Just 4.5") %>% bind_rows(endo.out6 %>% mutate(SHAPE="4.5 to 4.6"), endo.out7 %>% mutate(SHAPE="Just 4.6"))) %>%
ggplot()+ geom_line(aes(x=TAs, y=enbal.QGEN, group=SHAPE, colour=SHAPE), linewidth=2)+ theme_classic()
# Plotting all variables to look at patterns that may explain
for(i in 1:ncol(endo.out5)){
min5 <- min(endo.out5[, i], na.rm=TRUE)
max5 <- max(endo.out5[, i], na.rm=TRUE)
plot(TAs, endo.out5[, i], type="l", main=paste0("Variable ", colnames(endo.out5)[i]), ylim=c(min5*0.5, max5*1.5))
lines(TAs, endo.out6[, i], col="blue")
lines(TAs, endo.out7[, i], col="red")
abline(v=-4, col="grey")
abline(v=11, col="grey")
}
# From looking at the variables. The only one I can see that is a candidate is the K_FLESH
(endo.out5 %>% mutate(SHAPE="Just 4.5") %>% bind_rows(endo.out6 %>% mutate(SHAPE="4.5 to 4.6"), endo.out7 %>% mutate(SHAPE="Just 4.6"))) %>%
ggplot()+ geom_line(aes(x=TAs, y=treg.K_FLESH, group=SHAPE, colour=SHAPE), linewidth=2)+ theme_classic() +
geom_vline(xintercept=-4, colour="grey")+ geom_vline(xintercept=11, colour="grey")
# it appears to be caused by the delay in the increase in flesh conductivity
# Can see here that the transition in SHAPE_B occurs but the skin conductivity doesn't
endo.out6 %>%
ggplot()+ geom_line(aes(x=TAs, y=treg.SHAPE_B), linewidth=2)+ theme_classic() +
geom_vline(xintercept=-4, colour="grey")+ geom_vline(xintercept=11, colour="grey")
NicheMapR
Apr 16, 2026, 4:28:15 AMApr 16
to NicheMapR
Hi Shane,
Interesting issue - the calculation of the fat layer for an ellipsoid is quite a complex one. I tried the same code for a cylinder and the issue didn't happen. I'll dig into it soon and let you know.
All the best,
Mike
Mike
Shane Morris
Apr 21, 2026, 5:04:02 AMApr 21
to NicheMapR
Hi Mike,
Thanks, that would be great. Yes, it is much more complex than I initially thought. I was not aware that the fat thickness changed so much as the animal "stretches out".
Cheers,
Shane.
Shane Morris
7:47 AM (1 hour ago) 7:47 AM
to NicheMapR
Hi Mike,
Have you had a chance to look at this? Is there anyway to simply include a minimum fat thickness to avoid the fat from disappearing? As you can the fat thickness should likely be set to a minimum of 0.02 - 0.025.
Cheers,
Shane.
Reply all
Reply to author
Forward
0 new messages