Skip to content

Instantly share code, notes, and snippets.

@danlewer
Last active February 11, 2021 16:08
Show Gist options
  • Save danlewer/dcc13f0d01d2a0dd4c8266690927b9fa to your computer and use it in GitHub Desktop.
Save danlewer/dcc13f0d01d2a0dd4c8266690927b9fa to your computer and use it in GitHub Desktop.
# Calculate predicted Forced Exhaled Volume in 1 second (FEV1), based on age, sex, height, and ethnicity
# Function takes vectors
# Formula is from:
# ERS TASK FORCE REPORT. Multi-ethnic reference values for spirometry for the 3-95-yr age range: the global lung function 2012 equations
# Philip H.
# Source: Guideline 2012
# https://www.ers-education.org/lr/show-details/?idP=138497
fev1pred <- function (sex = 'f', heightCM = 180, ageYears = 50, eth = 'White') {
ind <- (sex == 'f') + 1
a0 <- c(-10.342, -9.6987)[ind] # intercept
a1 <- c(2.2196, 2.1211)[ind] # height in cm
a2 <- c(0.0574, -0.0270)[ind] # age
a3 <- c(-0.1589, -0.1484)[ind] # AfricanAmerican
a4 <- c(-0.0351, -0.0149)[ind] # NorthEastAsian
a5 <- c(-0.0881, -0.1208)[ind] # SouthEastAsian
a6 <- c(-0.0708, -0.0708)[ind] # OtherMixed
AfricanAmerican <- eth == 'AfricanAmerican'
NorthEastAsian <- eth == 'NorthEastAsian'
SouthEastAsian <- eth == 'SouthEastAsian'
OtherMixed <- eth == 'OtherMixed'
mspline_age <- c(3, 3.25, 3.5, 3.75, 4, 4.25, 4.5, 4.75, 5, 5.25, 5.5, 5.75, 6, 6.25, 6.5, 6.75, 7, 7.25, 7.5, 7.75, 8, 8.25, 8.5, 8.75, 9, 9.25, 9.5, 9.75, 10, 10.25, 10.5, 10.75, 11, 11.25, 11.5, 11.75, 12, 12.25, 12.5, 12.75, 13, 13.25, 13.5, 13.75, 14, 14.25, 14.5, 14.75, 15, 15.25, 15.5, 15.75, 16, 16.25, 16.5, 16.75, 17, 17.25, 17.5, 17.75, 18, 18.25, 18.5, 18.75, 19, 19.25, 19.5, 19.75, 20, 20.25, 20.5, 20.75, 21, 21.25, 21.5, 21.75, 22, 22.25, 22.5, 22.75, 23, 23.25, 23.5, 23.75, 24, 24.25, 24.5, 24.75, 25, 25.25, 25.5, 25.75, 26, 26.25, 26.5, 26.75, 27, 27.25, 27.5, 27.75, 28, 28.25, 28.5, 28.75, 29, 29.25, 29.5, 29.75, 30, 30.25, 30.5, 30.75, 31, 31.25, 31.5, 31.75, 32, 32.25, 32.5, 32.75, 33, 33.25, 33.5, 33.75, 34, 34.25, 34.5, 34.75, 35, 35.25, 35.5, 35.75, 36, 36.25, 36.5, 36.75, 37, 37.25, 37.5, 37.75, 38, 38.25, 38.5, 38.75, 39, 39.25, 39.5, 39.75, 40, 40.25, 40.5, 40.75, 41, 41.25, 41.5, 41.75, 42, 42.25, 42.5, 42.75, 43, 43.25, 43.5, 43.75, 44, 44.25, 44.5, 44.75, 45, 45.25, 45.5, 45.75, 46, 46.25, 46.5, 46.75, 47, 47.25, 47.5, 47.75, 48, 48.25, 48.5, 48.75, 49, 49.25, 49.5, 49.75, 50, 50.25, 50.5, 50.75, 51, 51.25, 51.5, 51.75, 52, 52.25, 52.5, 52.75, 53, 53.25, 53.5, 53.75, 54, 54.25, 54.5, 54.75, 55, 55.25, 55.5, 55.75, 56, 56.25, 56.5, 56.75, 57, 57.25, 57.5, 57.75, 58, 58.25, 58.5, 58.75, 59, 59.25, 59.5, 59.75, 60, 60.25, 60.5, 60.75, 61, 61.25, 61.5, 61.75, 62, 62.25, 62.5, 62.75, 63, 63.25, 63.5, 63.75, 64, 64.25, 64.5, 64.75, 65, 65.25, 65.5, 65.75, 66, 66.25, 66.5, 66.75, 67, 67.25, 67.5, 67.75, 68, 68.25, 68.5, 68.75, 69, 69.25, 69.5, 69.75, 70, 70.25, 70.5, 70.75, 71, 71.25, 71.5, 71.75, 72, 72.25, 72.5, 72.75, 73, 73.25, 73.5, 73.75, 74, 74.25, 74.5, 74.75, 75, 75.25, 75.5, 75.75, 76, 76.25, 76.5, 76.75, 77, 77.25, 77.5, 77.75, 78, 78.25, 78.5, 78.75, 79, 79.25, 79.5, 79.75, 80, 80.25, 80.5, 80.75, 81, 81.25, 81.5, 81.75, 82, 82.25, 82.5, 82.75, 83, 83.25, 83.5, 83.75, 84, 84.25, 84.5, 84.75, 85, 85.25, 85.5, 85.75, 86, 86.25, 86.5, 86.75, 87, 87.25, 87.5, 87.75, 88, 88.25, 88.5, 88.75, 89, 89.25, 89.5, 89.75, 90, 90.25, 90.5, 90.75, 91, 91.25, 91.5, 91.75, 92, 92.25, 92.5, 92.75, 93, 93.25, 93.5, 93.75, 94, 94.25, 94.5, 94.75, 95)
mspline <- cbind(male = c(-0.1133,-0.1073,-0.1011,-0.0951,-0.0893,-0.0841,-0.0799,-0.0769,-0.0752,-0.075,-0.0758,-0.0771,-0.0787,-0.0803,-0.0816,-0.0823,-0.0822,-0.0815,-0.0804,-0.0792,-0.0778,-0.0763,-0.0745,-0.0721,-0.0691,-0.0658,-0.0622,-0.0586,-0.0549,-0.0513,-0.0476,-0.0437,-0.0395,-0.035,-0.0299,-0.0241,-0.0176,-0.0101,-0.0019,0.0071,0.0169,0.0274,0.0384,0.0497,0.0612,0.0728,0.0844,0.0958,0.1068,0.1175,0.1276,0.1371,0.146,0.1542,0.1616,0.1684,0.1744,0.1798,0.1845,0.1887,0.1924,0.1956,0.1984,0.2008,0.2029,0.2046,0.206,0.2072,0.2081,0.2087,0.209,0.2092,0.2091,0.2089,0.2084,0.2079,0.2071,0.2063,0.2053,0.2042,0.203,0.2016,0.2002,0.1987,0.197,0.1954,0.1936,0.1918,0.1899,0.188,0.1861,0.1841,0.1821,0.1801,0.1781,0.176,0.1739,0.1718,0.1697,0.1677,0.1656,0.1635,0.1615,0.1594,0.1574,0.1554,0.1534,0.1514,0.1495,0.1475,0.1455,0.1436,0.1417,0.1397,0.1378,0.1359,0.134,0.1321,0.1302,0.1283,0.1265,0.1246,0.1227,0.1209,0.119,0.1172,0.1153,0.1135,0.1116,0.1097,0.1078,0.1059,0.104,0.1021,0.1001,0.0982,0.0962,0.0943,0.0923,0.0903,0.0883,0.0863,0.0843,0.0823,0.0803,0.0782,0.0762,0.0742,0.0721,0.07,0.068,0.0659,0.0638,0.0617,0.0596,0.0575,0.0554,0.0533,0.0511,0.049,0.0469,0.0448,0.0427,0.0406,0.0386,0.0365,0.0344,0.0323,0.0302,0.0281,0.0261,0.024,0.0219,0.0198,0.0177,0.0156,0.0135,0.0114,0.0093,0.0072,0.005,0.0029,0.0007,-0.0015,-0.0036,-0.0058,-0.008,-0.0103,-0.0125,-0.0147,-0.017,-0.0193,-0.0216,-0.0239,-0.0262,-0.0285,-0.0309,-0.0332,-0.0356,-0.038,-0.0404,-0.0428,-0.0453,-0.0478,-0.0502,-0.0527,-0.0552,-0.0578,-0.0603,-0.0629,-0.0654,-0.068,-0.0706,-0.0732,-0.0759,-0.0785,-0.0812,-0.0839,-0.0866,-0.0893,-0.092,-0.0947,-0.0975,-0.1002,-0.103,-0.1058,-0.1086,-0.1114,-0.1143,-0.1171,-0.1199,-0.1228,-0.1257,-0.1286,-0.1315,-0.1344,-0.1373,-0.1402,-0.1431,-0.1461,-0.149,-0.1519,-0.1549,-0.1578,-0.1608,-0.1638,-0.1667,-0.1697,-0.1727,-0.1757,-0.1786,-0.1816,-0.1846,-0.1876,-0.1906,-0.1936,-0.1966,-0.1996,-0.2026,-0.2056,-0.2086,-0.2116,-0.2147,-0.2177,-0.2207,-0.2237,-0.2267,-0.2298,-0.2328,-0.2358,-0.2388,-0.2418,-0.2449,-0.2479,-0.2509,-0.2539,-0.2569,-0.2599,-0.263,-0.266,-0.269,-0.272,-0.275,-0.278,-0.281,-0.284,-0.2869,-0.2899,-0.2929,-0.2959,-0.2989,-0.3018,-0.3048,-0.3077,-0.3107,-0.3136,-0.3166,-0.3195,-0.3224,-0.3253,-0.3282,-0.3311,-0.334,-0.3369,-0.3398,-0.3427,-0.3455,-0.3484,-0.3512,-0.3541,-0.3569,-0.3597,-0.3625,-0.3654,-0.3682,-0.3709,-0.3737,-0.3765,-0.3793,-0.382,-0.3848,-0.3875,-0.3903,-0.393,-0.3957,-0.3984,-0.4011,-0.4038,-0.4065,-0.4092,-0.4119,-0.4145,-0.4172,-0.4198,-0.4225,-0.4251,-0.4277,-0.4303,-0.4329,-0.4355,-0.4381,-0.4407,-0.4433,-0.4459,-0.4484,-0.451,-0.4536,-0.4561,-0.4586,-0.4612,-0.4637,-0.4662,-0.4687,-0.4712,-0.4737,-0.4762,-0.4787,-0.4811,-0.4836,-0.4861,-0.4885,-0.491,-0.4934,-0.4959,-0.4983,-0.5007,-0.5031,-0.5055,-0.5079),
female = c(-0.2311,-0.217,-0.204,-0.1922,-0.1817,-0.1727,-0.1651,-0.1592,-0.1548,-0.1518,-0.1494,-0.1474,-0.1452,-0.1426,-0.1393,-0.1354,-0.131,-0.1264,-0.1217,-0.1171,-0.1125,-0.1076,-0.1022,-0.0963,-0.0897,-0.0825,-0.0747,-0.0663,-0.0573,-0.0479,-0.038,-0.0278,-0.0172,-0.0063,0.0048,0.0161,0.0274,0.0386,0.0496,0.0604,0.0709,0.081,0.0907,0.0999,0.1086,0.1168,0.1244,0.1315,0.1379,0.1438,0.1492,0.154,0.1583,0.1621,0.1655,0.1684,0.1711,0.1733,0.1753,0.177,0.1785,0.1797,0.1808,0.1816,0.1823,0.1829,0.1833,0.1837,0.1839,0.1841,0.1842,0.1842,0.1841,0.184,0.1838,0.1835,0.1832,0.1828,0.1823,0.1818,0.1812,0.1806,0.1799,0.1792,0.1785,0.1777,0.1769,0.1761,0.1753,0.1745,0.1737,0.1729,0.1721,0.1713,0.1705,0.1697,0.169,0.1682,0.1674,0.1666,0.1658,0.165,0.1642,0.1634,0.1625,0.1617,0.1608,0.1599,0.159,0.1581,0.1572,0.1562,0.1553,0.1543,0.1533,0.1523,0.1512,0.1501,0.149,0.1479,0.1467,0.1456,0.1444,0.1431,0.1418,0.1406,0.1392,0.1379,0.1365,0.1351,0.1337,0.1322,0.1308,0.1292,0.1277,0.1262,0.1246,0.123,0.1214,0.1197,0.118,0.1164,0.1147,0.1129,0.1112,0.1094,0.1076,0.1058,0.104,0.1022,0.1003,0.0985,0.0966,0.0947,0.0928,0.0909,0.0889,0.087,0.085,0.083,0.0811,0.0791,0.0771,0.0751,0.0731,0.071,0.069,0.067,0.065,0.063,0.0609,0.0589,0.0568,0.0548,0.0527,0.0507,0.0486,0.0465,0.0445,0.0424,0.0403,0.0382,0.0361,0.0339,0.0318,0.0297,0.0275,0.0254,0.0232,0.021,0.0188,0.0166,0.0144,0.0122,0.0099,0.0077,0.0054,0.0032,0.0009,-0.0014,-0.0037,-0.0061,-0.0084,-0.0108,-0.0131,-0.0155,-0.0179,-0.0203,-0.0227,-0.0252,-0.0276,-0.0301,-0.0326,-0.035,-0.0375,-0.0401,-0.0426,-0.0451,-0.0477,-0.0503,-0.0529,-0.0555,-0.0581,-0.0607,-0.0634,-0.066,-0.0687,-0.0714,-0.0741,-0.0768,-0.0795,-0.0822,-0.085,-0.0878,-0.0905,-0.0933,-0.0961,-0.0989,-0.1018,-0.1046,-0.1075,-0.1103,-0.1132,-0.1161,-0.119,-0.1219,-0.1249,-0.1278,-0.1308,-0.1338,-0.1368,-0.1398,-0.1428,-0.1458,-0.1488,-0.1519,-0.155,-0.158,-0.1611,-0.1642,-0.1674,-0.1705,-0.1736,-0.1768,-0.1799,-0.1831,-0.1863,-0.1895,-0.1926,-0.1958,-0.1991,-0.2023,-0.2055,-0.2087,-0.212,-0.2152,-0.2184,-0.2217,-0.2249,-0.2282,-0.2315,-0.2347,-0.238,-0.2413,-0.2445,-0.2478,-0.2511,-0.2543,-0.2576,-0.2609,-0.2642,-0.2674,-0.2707,-0.274,-0.2773,-0.2805,-0.2838,-0.2871,-0.2903,-0.2936,-0.2968,-0.3001,-0.3033,-0.3065,-0.3098,-0.313,-0.3162,-0.3194,-0.3226,-0.3258,-0.329,-0.3322,-0.3354,-0.3386,-0.3417,-0.3449,-0.348,-0.3512,-0.3543,-0.3574,-0.3606,-0.3637,-0.3668,-0.3699,-0.373,-0.376,-0.3791,-0.3822,-0.3852,-0.3883,-0.3913,-0.3944,-0.3974,-0.4004,-0.4034,-0.4064,-0.4094,-0.4124,-0.4153,-0.4183,-0.4213,-0.4242,-0.4272,-0.4301,-0.433,-0.4359,-0.4389,-0.4418,-0.4446,-0.4475,-0.4504,-0.4533,-0.4561,-0.459,-0.4618,-0.4647,-0.4675,-0.4703,-0.4732,-0.476,-0.4788,-0.4816,-0.4844,-0.4871,-0.4899,-0.4927,-0.4954,-0.4982,-0.5009))
M <- mspline[cbind(findInterval(ageYears, mspline_age), ind)]
exp(a0 + a1*log(heightCM) + a2*log(ageYears) + a3*AfricanAmerican + a4*NorthEastAsian + a5*SouthEastAsian + a6*OtherMixed + M)
}
# examples in guidance
fev1pred(sex = c('m', 'm', 'm', 'f', 'f'),
heightCM = c(107, 152, 175, 165, 165),
ageYears = c(4.8, 12.2, 53, 39.1, 39.1),
eth = c('white', 'AfricanAmerican', 'SouthEastAsian', 'SouthEastAsian', 'OtherMixed'))
# given values in guidance: 1.04, 2.19, 3.39, 2.78, 2.92
# check efficiency
reps <- 1e6
system.time(
fev1pred(sex = sample(c('m', 'f'), reps, T),
heightCM = sample(130:210, reps, T),
ageYears = rnorm(reps, 50, 7),
eth = sample(c('white', 'AfricanAmerican', 'SouthEastAsian', 'SouthEastAsian', 'OtherMixed'), reps, T))
)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment