Estimates by samples tables hard to interpret

[erex]

Prompted by this question to the list 04Oct22, I replicate the problem with data set LTExercise distributed with the package. Essence of problem: covered area for each transect is computed as constant even though transect lengths differ.

I looked at one of the "tame" data sets that ships with the Distance package, LTExercise​: line transect data with 12 replicate transects, true density 79.8 animals per km^2, just to have a look at some results.

library(Distance)
data(LTExercise)
conv <- convert_units("meter", "kilometer", "square kilometer")
lt_hn <- ds(LTExercise, key="hn", convert_units = conv)
summary(lt_hn)

                       Estimate          SE         CV
Average p             0.3491863  0.02160948 0.06188525
N in covered region 300.6991117 30.11199438 0.10013995

Summary statistics:
   Region Area CoveredArea Effort   n  k     ER     se.ER     cv.ER
1 Default    1      3.4368     48 105 12 2.1875 0.3169604 0.1448962

Density estimate for the study area:

Density:
  Label Estimate       se        cv      lcl      ucl       df
1 Total 87.49392 13.78541 0.1575585 62.70329 122.0859 15.32522

Density estimates by transect

> lt_hn$dht$individuals$bysample
    Region  Sample Effort Sample.Area Area  n      Nhat     Nchat      Dhat
1  Default  Line 1      5      3.4368    1 14 11.665857 40.093215 11.665857
2  Default Line 10      7      3.4368    1 23 19.165336 65.867424 19.165336
3  Default Line 11      3      3.4368    1  0  0.000000  0.000000  0.000000
4  Default Line 12      4      3.4368    1  5  4.166377 14.319005  4.166377
5  Default  Line 2      2      3.4368    1  3  2.499826  8.591403  2.499826
6  Default  Line 3      6      3.4368    1  8  6.666204 22.910409  6.666204
7  Default  Line 4      4      3.4368    1  9  7.499479 25.774210  7.499479
8  Default  Line 5      3      3.4368    1  3  2.499826  8.591403  2.499826
9  Default  Line 6      1      3.4368    1  4  3.333102 11.455204  3.333102
10 Default  Line 7      4      3.4368    1 11  9.166030 31.501812  9.166030
11 Default  Line 8      4      3.4368    1  8  6.666204 22.910409  6.666204
12 Default  Line 9      5      3.4368    1 17 14.165683 48.684618 14.165683

Note that for the estimates by transect, the Sample.Area is reported to be the same, even though the transects differ in length. I think this might be the source of the discrepancy you describe.

[LHMarshall]

This is not a bug but it is a table that could do with improvement! We should add this to discussion at the next development meeting @lenthomas @erex We need to decide what is useful information to the user, at the moment it is unclear how to interpret the provided output.

This table currently provides the following information:

• Region: region / stratum label
• Sample: sample label
• Effort: line length of the sampler
• Sample.Area - Covered Area in the region / stratum (this seems a misleading name however like it should be the sampled area for the sampler)
• n: number of observations for the sampler
• Nhat: abundance contribution for the observations on the sampler to the overall study region Nchat*(A/Ac)
  [summing the Nhat column gives the total abundance estimate]
• Nchat: abundance contribution for the observations on the sampler to the overall covered area (i.e. n_i/Pa)
  [i.e. Nhat = Nchat*Area/Covered Area for each stratum]
• Dhat: Nchat/Covered Area
  [summing Dhat gives overall density estimates]

See code below for calculation verification (note I have changes the overall area to be larger than the covered region and it also makes the calculations more obvious as having an area of 1 was a special case)

library(Distance)
data(LTExercise)
dht.dat <- unflatten(LTExercise)
conv <- convert_units("meter", "kilometer", "square kilometer")

lt_ddf <- ddf(dsmodel=~cds(key = "hn", formula = ~1),
              data = LTExercise)

dht.dat$region.table$Area <- 5.75

lt_dht <- dht(lt_ddf, 
              region.table = dht.dat$region.table,
              sample.table = dht.dat$sample.table, 
              obs.table = dht.dat$obs.table,
              options = list(convert.units = conv))

> print(lt_dht, bysample = TRUE)
Abundance and density estimates from distance sampling
Variance       : R2, N/L 

Summary statistics

   Region Area CoveredArea Effort   n  k     ER     se.ER     cv.ER
1 Default 5.75      3.4368     48 105 12 2.1875 0.3169604 0.1448962

Abundance:
  Region Estimate       se        cv      lcl      ucl       df
1  Total 503.0901 79.26612 0.1575585 360.5439 701.9939 15.32522

Density:
  Region Estimate       se        cv      lcl      ucl       df
1  Total 87.49392 13.78541 0.1575585 62.70329 122.0859 15.32522

Estimates by sample:
    Region  Sample Effort Sample.Area Area  n      Nhat     Nchat      Dhat
1  Default  Line 1      5      3.4368 5.75 14  67.07868 40.093215 11.665857
2  Default Line 10      7      3.4368 5.75 23 110.20068 65.867424 19.165336
3  Default Line 11      3      3.4368 5.75  0   0.00000  0.000000  0.000000
4  Default Line 12      4      3.4368 5.75  5  23.95667 14.319005  4.166377
5  Default  Line 2      2      3.4368 5.75  3  14.37400  8.591403  2.499826
6  Default  Line 3      6      3.4368 5.75  8  38.33067 22.910409  6.666204
7  Default  Line 4      4      3.4368 5.75  9  43.12201 25.774210  7.499479
8  Default  Line 5      3      3.4368 5.75  3  14.37400  8.591403  2.499826
9  Default  Line 6      1      3.4368 5.75  4  19.16534 11.455204  3.333102
10 Default  Line 7      4      3.4368 5.75 11  52.70467 31.501812  9.166030
11 Default  Line 8      4      3.4368 5.75  8  38.33067 22.910409  6.666204
12 Default  Line 9      5      3.4368 5.75 17  81.45268 48.684618 14.165683
# Use Dhats to calculate overall density
> sum(lt_dht$individuals$bysample$Dhat)
[1] 87.49392
# Use Nhats to calculate overall abundance
> sum(lt_dht$individuals$bysample$Nhat)
[1] 503.0901
# Calculate Nchat
> (lt_dht$individuals$bysample$n/0.3491863)
 [1] 40.093211 65.867418  0.000000 14.319004  8.591402 22.910406 25.774207  8.591402 11.455203 31.501809 22.910406 48.684613
# Calculate Nhat
> lt_dht$individuals$bysample$Nchat*(lt_dht$individuals$bysample$Area/lt_dht$individuals$bysample$Sample.Area)
 [1]  67.07868 110.20068   0.00000  23.95667  14.37400  38.33067  43.12201  14.37400  19.16534  52.70467  38.33067  81.45268

[erex]

If Sample.Area and Area are not transect specific, perhaps they should not be included in the table.

[LHMarshall]

Create a new more user friendly table for output with:
Sample.Area - make this now the covered area of the transect
Nchat - N in covered area of the transect
Dchat - D in covered area of the transect