Spatio-temporal analysis of H5N1 in cats in Poland (comparision of sources)

2023-08-05

Author: Andrzej Jarynowski, IVEB FU Berlin

1. Load the data and introduction

Loading and preparing the unofficial data - suspected 87 case submitted by animal owners. Loading and preparing the official data (from reference lab) in an area of Poland. It consists of 30 cases (positive results) and 27 controls (negative results).

1.1 Can we trust unofficial sources?

The total population of on farms (~2M) and domestic cats (~2M) in Poland is about 4 million (https://europeanpetfood.org/wp-content/uploads/2023/02/Annual-Report-2022-2.pdf , https://www.sciencedirect.com/science/article/pii/S2351989418303196?via%3Dihub) which shows that less than every third household in Poland has a cat.

Thus means that the density of cats per number of inhabitants is about 1.5-2 times lower than the European average (https://landgeist.com/2023/04/29/cat-ownership-in-europe/).

1.2 Media listening

We added geographical distribution of queries for cats health and AI related topics. Digital epidemiology approach is widely used, for example, “infoveillance” as a term comes from information surveillance. Variable “interest” is mean across AI, cats disease and cats deaths queries.

2. Comparison

Comparative analysis and interpretation of data to determine:

  • Triangulated distribution of diseases in time and space

  • Presence or absence of disease

  • source detection

2.1 Spatial component

par(mfrow=c(1,3))

colourScale = colorRampPalette(brewer.pal(9,"OrRd"))(16)[1:16]
cols = colourScale[zlicz$cases_per_million+1]
plot(contour, lwd=0.4, border="gray30", col=NA); plot(nutsM, col=cols , lwd=0.1, add=T)
rast = raster(as.matrix(c(min(zlicz$cases_per_million,na.rm=T),max(zlicz$cases_per_million,na.rm=T))))
title("Suspected cases per mil cit", side=3, line=-1.5, at=0, cex=0.50, col="gray30")
## Warning in title("Suspected cases per mil cit", side = 3, line = -1.5, at = 0,
## : 'side' nie jest parametrem graficznym
## Warning in title("Suspected cases per mil cit", side = 3, line = -1.5, at = 0,
## : 'at' nie jest parametrem graficznym
plot(rast, legend.only=T, add=T, col=colourScale, legend.width=0.5, legend.shrink=0.3, smallplot=c(0.870,0.885,0.3,0.8), adj=3,
     axis.args=list(cex.axis=0.65, lwd=0, col="gray30", lwd.tick=0.2, col.tick="gray30", tck=-1.0, col.axis="gray30", line=0, mgp=c(0,0.45,0)), alpha=1, side=3)
colourScale = colorRampPalette(brewer.pal(9,"OrRd"))(16)[1:16]
cols = colourScale[woj_all_media$cases_per_million_r+1]
plot(contour, lwd=0.4, border="gray30", col=NA); plot(nutsM, col=cols , lwd=0.1, add=T)
rast = raster(as.matrix(c(min(woj_all_media$cases_per_million_r,na.rm=T),max(woj_all_media$cases_per_million_r,na.rm=T))))
title("Confirmed cases per mil cit",   title.cex = 0.2, side=3, line=-1.5, at=0, cex=0.50, col="gray30")
## Warning in title("Confirmed cases per mil cit", title.cex = 0.2, side = 3, :
## 'title.cex' nie jest parametrem graficznym
## Warning in title("Confirmed cases per mil cit", title.cex = 0.2, side = 3, :
## 'side' nie jest parametrem graficznym
## Warning in title("Confirmed cases per mil cit", title.cex = 0.2, side = 3, :
## 'at' nie jest parametrem graficznym
plot(rast, legend.only=T, add=T, col=colourScale, legend.width=0.5, legend.shrink=0.3, smallplot=c(0.870,0.885,0.3,0.8), adj=3,
     axis.args=list(cex.axis=0.65, lwd=0, col="gray30", lwd.tick=0.2, col.tick="gray30", tck=-1.0, col.axis="gray30", line=0, mgp=c(0,0.45,0)), alpha=1, side=3)
colourScale = colorRampPalette(brewer.pal(9,"OrRd"))(16)[1:16]
cols = colourScale[woj_all_media$interest/15+1]
plot(contour, lwd=0.4, border="gray30", col=NA); plot(nutsM, col=cols , lwd=0.1, add=T)
rast = raster(as.matrix(c(min(woj_all_media$interest/15,na.rm=T),max(woj_all_media$interest/15,na.rm=T))))
title("Interest in Google",   title.cex = 0.2, side=3, line=-1.5, at=0, cex=0.50, col="gray30")
## Warning in title("Interest in Google", title.cex = 0.2, side = 3, line = -1.5,
## : 'title.cex' nie jest parametrem graficznym
## Warning in title("Interest in Google", title.cex = 0.2, side = 3, line = -1.5,
## : 'side' nie jest parametrem graficznym
## Warning in title("Interest in Google", title.cex = 0.2, side = 3, line = -1.5,
## : 'at' nie jest parametrem graficznym
plot(rast, legend.only=T, add=T, col=colourScale, legend.width=0.5, legend.shrink=0.3, smallplot=c(0.870,0.885,0.3,0.8), adj=3,
     axis.args=list(cex.axis=0.65, lwd=0, col="gray30", lwd.tick=0.2, col.tick="gray30", tck=-1.0, col.axis="gray30", line=0, mgp=c(0,0.45,0)), alpha=1, side=3)

We see clear overall between all datasets in Pomerania region, which definitely was the most affected.

newmap <- getMap(resolution = "low")
plot(newmap, xlim = c(15, 25), ylim = c(45, 55), asp = 1)
points(koty_sorter$lon+rnorm(57)/10, koty_sorter$lat+rnorm(57)/10, col = koty_sorter$case+1, pch=20, cex = 1)
points(cats_geo$lon +rnorm(78)/10, cats_geo$lat +rnorm(78)/10, col = 3, pch=20, cex = 1, add=T)
## Warning in plot.xy(xy.coords(x, y), type = type, ...): 'add' nie jest
## parametrem graficznym
legend("topleft", legend = c("-","+","citizen"), col = c(1:3), lty = 1, bty = "n")

title("Cases (+) and control (-) and participatory epidemiology ")

2.2 Temporal component

weekly_plus=ddply(koty_plus,.(week), summarise, cases_plus=length(lat))
weekly_sus=ddply(cats,.(week), summarise, cases_citizen=length(Województwo))
weekly_sus$cases_plus=NA
weekly_sus[9,1]=28
weekly_sus$cases_plus[9]=5
weekly_sus$cases_plus[8]=4
weekly_sus$cases_plus[7]=9
weekly_sus$cases_plus[6]=12
weekly_gt=ddply(gt,.(week), summarise, case_gt=mean(interest))
weekly_sus2=left_join(weekly_sus,weekly_gt,by=c("week"="week"))
media_interest$date=as.Date(media_interest$Day)
media_interest$week=week(media_interest$date)
weekly_media=ddply(media_interest,.(week), summarise, case_media=mean(Mentions.count))
weekly_sus2=left_join(weekly_sus2,weekly_media,by=c("week"="week"))


ccfvalues = ccf(weekly_sus$cases_plus,weekly_sus$cases_citizen, na.action = na.pass, lag.max=2)

ccfvalues = ccf(weekly_sus2$cases_plus,weekly_sus2$case_gt, na.action = na.pass, lag.max=2)

ccfvalues = ccf(weekly_sus2$cases_plus,weekly_sus2$case_media, na.action = na.pass, lag.max=2)

Citizens reports are one week ahead the positive cases. Interest in Google trends and in media is one week after positive cases.

daily_plus=ddply(koty_plus,.(date), summarise, cases_plus=length(lat))
daily_sus=ddply(cats,.(date), summarise, cases_citizen=length(Województwo))
 inter=inter[,c(1:3)]
 inter$date=as.Date(inter$date)
inter$media_interest=0
inter$google_trends=0
inter$media_interest[12:40]=media_interest$Mentions.count[1:29]/390
inter$google_trends[12:40]=gt$interest[1:29]/10

por_long <- melt(inter, id="date") 
por_long$date=as.Date(por_long$date)
Sys.setlocale(category="LC_ALL", locale="US")
## Warning in Sys.setlocale(category = "LC_ALL", locale = "US"): using locale code
## page other than 65001 ("UTF-8") may cause problems
## [1] "LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252"
p = ggplot(por_long ,aes(date,value))+ geom_line(aes(color = variable, linetype = variable)) +ggtitle("Interpolated daily series of signals")+   theme_bw()
p

Sys.setlocale(category="LC_ALL", locale="Poland")
## Warning in Sys.setlocale(category = "LC_ALL", locale = "Poland"): Zadania
## raportów OS aby ustawic lokalizacje na "Poland" nie moga zostac wykonane
## [1] ""

Interpolated daily dynamics of positive case, suspected cases from participatory study, normalized interest in Google trends and No articles/post on the issue per million citizen in media.

knitr::include_graphics("przebieg_events.png", dpi=100)

Trajectory of main events in outbreak investigation. Based on all sources (all of which have their advantages and disadvantages), it seems that the epizootic has already passed the peak phase and is in the ending phase.

kor_=cor(as.matrix(inter[,-1]), use="pairwise.complete.obs")

pheatmap(kor_, display_numbers = T)

Correlations between time series showed that there is de-synchronization between participatory reporting and supply/demand information needs.

for (i in 2:4)
{ for (j in (i+1):5)
{
  

  
ccfvalues = ccf(inter[,i], inter[,j], na.action = na.pass, lag.max=21)
print(colnames(inter)[i])
print(colnames(inter)[j])
}}

## [1] "cases_citizen"
## [1] "cases_plus"

## [1] "cases_citizen"
## [1] "media_interest"

## [1] "cases_citizen"
## [1] "google_trends"

## [1] "cases_plus"
## [1] "media_interest"

## [1] "cases_plus"
## [1] "google_trends"

## [1] "media_interest"
## [1] "google_trends"

Lag-Correlations between time series in days. Participatory reporting is aheading supply/demand information needs. There is interesting bimodal strong lagged correlation between positive cases and No. of mentions in media. There is one peak of correlation, that media report information one day after testing positive cases (mainly local media ones based on information from local veterinary authorities) and one week before.

There is another bimodal medium lagged correlation between positive cases and participatory reporting. There is one peak of correlation one week and one 2 before.

Sys.setlocale(category="LC_ALL", locale="US")
## Warning in Sys.setlocale(category = "LC_ALL", locale = "US"): using locale code
## page other than 65001 ("UTF-8") may cause problems
## [1] "LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252"
cats$cum=NA
cats = cats%>% mutate(cum = cumsum(name))
cats2=cats
cats2[88,28]=koty_plus$date[30]
cats2$cum[88]=cats2$cum[87]
ggplot()+geom_line(data=cats2,aes(date,cum), color="blue")+geom_line(data=koty_plus,aes(date,cum),color="black")+
xlab("Time")+ylab("Cumulative suspected cases")+
  theme_bw()

koty_plus2=koty_plus
koty_plus2$code="plus"
koty_minus=koty_sorter[koty_sorter$case==0,]
koty_minus$code="minus"
cats2=cats
cats2$code="citizens"
dat <- rbind(cats2[,c(28,34)],koty_plus2[,c(9,18)],koty_minus[,c(9,18)])
por_long <- melt(dat, id="code") 
px<-ggplot(por_long, aes(x=as.factor(code), y=value, fill=as.factor(code))) +
  geom_violin(trim=FALSE, show.legend = TRUE)+ ggtitle("Smoothed positive cases and suspeceted (from citizens) intensity in time")+ ylab("time")+
  theme_bw()
px

## ANN citizens
mean(nndist(pp_, k = 1))
## [1] 0.2038613
mean(nndist(pp_, k = 2))
## [1] 0.3584615
## ANN laboratory plus

mean(nndist(cat_case, k = 1))
## [1] 0.1825677
mean(nndist(cat_case, k = 2))
## [1] 0.3020286
## ANN laboratory minus
mean(nndist(cat_cont, k = 1))
## [1] 0.3774133
mean(nndist(cat_cont, k = 2))
## [1] 0.6496076

Positive data points has smaller ANN than suspected (even there is much less of them), thus is suggest that positive cases form some characteristics chains. This could support migratory bird cause of outbreak.

2.3 Spatio-Temporal correlation

Let look at hierarchical correlation structure of series of provinces (where positive cases ending with r, and suspected from participatory sample are pure names). Mean/Medium Day is central date in the given province.

mat_char=woj_all_media[,c(3,4,9:11,13,17)]

mat_num = as.data.frame(sapply(mat_char, as.numeric))

kor2=cor(as.matrix(mat_num), use="pairwise.complete.obs")


pheatmap(kor2, display_numbers = T)

We see that interest in media are strictly correlated with participatory engagement.

Moreover, the most common geographical places till so called risk assessment phase mentioned in the media were Gdańsk, Lublin, Poznań, Bygdoszcz.

plus=which(koty_sorter$case==1)

par(mfrow=c(1,2))
newmap <- getMap(resolution = "low")

km.res <- kmeans(norm_scale[plus, c(1,2,3,3)], 3, nstart = 5)
#km.res$cluster
plot(newmap, xlim = c(15, 25), ylim = c(45, 55), asp = 1)
points(koty_plus$lon +rnorm(30)/10, koty_plus$lat +rnorm(30)/10, col = km.res$cluster, pch=20, cex = 1)
title("XYTT reference lab +")

newmap <- getMap(resolution = "low")

km.res <- kmeans(norm_scale_[, c(1,2,3,3)], 3, nstart = 5)
#km.res$cluster
plot(newmap, xlim = c(15, 25), ylim = c(45, 55), asp = 1)
points(cats_geo$lon +rnorm(78)/10, cats_geo$lat +rnorm(78)/10, col = km.res$cluster, pch=20, cex = 1)
title("XYTT participatory sample")

cats_geo$clus=km.res$cluster

summ_clus=ddply(cats_geo,.(clus), summarise, day_med=median(date), day_mean=mean(date),count=length(date), lon=mean(lon), lat=mean(lat))
summ_clus
##   clus    day_med   day_mean count      lon      lat
## 1    1 2023-06-17 2023-06-18    44 18.77591 51.60193
## 2    2 2023-05-20 2023-05-21     8 21.77408 50.99285
## 3    3 2023-06-17 2023-06-16    26 18.64486 53.99787

There is clear spatio-temporal cluster of suspected cases in second part of Mai on Lublin/Podkarpacie area. Time gives no additional information (in terms of quality of clustering) in reference lab positive cases.

2.4 Genetic distances

19 samples from WOAH reference lab proved that there are at least 2 separate sources of infections (one in Western Poland and second dispersed around the country) https://www.eurosurveillance.org/content/10.2807/1560-7917.ES.2023.28.31.2300366

l=str_remove_all(colnames(snp_1),"A/domestic cat/Poland/")
pheatmap(p_piwet, display_numbers = T,labels_row = l)

M1 <- mst(snp_1)
colnames(snp_1)
##  [1] "A/domestic cat/Poland/H246-M/2023 H5N1 2023-06-21 Gdynia brain" 
##  [2] "A/domestic cat/Poland/H255-M/2023 H5N1 2023-06-21 Pruszcz brain"
##  [3] "A/domestic cat/Poland/H256-G/2023 H5N1 2023-06-24 Lublin"       
##  [4] "A/domestic cat/Poland/H257-G/2023 H5N1 2023-06-24 Lublin"       
##  [5] "A/domestic cat/Poland/H263-G/2023 H5N1 2023-06-26 Komarow-Osada"
##  [6] "A/domestic cat/Poland/H264-G/2023 H5N1 2023-06-26 Poznan"       
##  [7] "A/domestic cat/Poland/H266-W/2023 H5N1 2023-06-19 Bydgoszcz"    
##  [8] "A/domestic cat/Poland/H267-W/2023 H5N1 Strzelin"                
##  [9] "A/domestic cat/Poland/H270-W/2023 H5N1 lubelskie"               
## [10] "A/domestic cat/Poland/H271-W/2023 H5N1 lubelskie"               
## [11] "A/domestic cat/Poland/H277-W1/2023 H5N1 2023-06-26 Namyslow"    
## [12] "A/domestic cat/Poland/H249/2023 H5N1 2023-06-22 Gdansk"         
## [13] "A/domestic cat/Poland/H248/2023 H5N1 2023-06-15 Pruszcz Gdanski"
## [14] "A/domestic cat/Poland/Kot2/2023 H5N1 2023-06-19 Poznan"         
## [15] "A/domestic cat/Poland/Kot1/2023 H5N1 2023-06-19 Poznan"         
## [16] "A/domestic cat/Poland/H254/2023 H5N1 2023-06-22 Lublin"         
## [17] "A/domestic cat/Poland/H252/2023 H5N1 2023-06-22 Lublin"         
## [18] "A/domestic cat/Poland/H247/2023 H5N1 2023-06-20 Gdansk"         
## [19] "A/domestic cat/Poland/H253/2023 H5N1 2023-06-22 Lublin"
plot(M1)

Node number is representing consequential samples in munimum spanning tree representation.

11 samples from https://www.eurosurveillance.org/content/10.2807/1560-7917.ES.2023.28.31.2300390#html_fulltext

pheatmap(p_pyrc, display_numbers = T,labels_row = colnames(snp_2))

M2 <- mst(snp_2)
colnames(snp_2)
##  [1] "A.environment.Poland.Kra1.2023"    "A.domestic_cat.Poland.H248.2023"  
##  [3] "A.domestic_cat.Poland.H252.2023"   "A.domestic_cat.Poland.H254.2023"  
##  [5] "A.domestic_cat.Poland.H247.2023"   "A.domestic_cat.Poland.H246.M.2023"
##  [7] "A.domestic_cat.Poland.H255.M.2023" "A.domestic_cat.Poland.H253.2023"  
##  [9] "A.domestic_cat.Poland.H249.2023"   "A.cat.Poland.Gda1.2023"           
## [11] "A.cat.Poland.Kra1.2023."
plot(M2, graph = "nsca")

Node number is representing consequential samples in munimum spanning tree representation.

Samples 1 [1] “A.environment.Poland.Kra1.2023” and 11 “A.cat.Poland.Kra1.2023.” are from suspected cat case and its feed (poultry meat), which never been confimed by WOAH reference lab.

As both articles with sequences have already appeared in Eurosurveillance, you can add a genetic dimension to the analysis. I would not agree with the authors that there are no spatial clusters. There are clearly two separate clusters for Wielkopolska and Pomerania. The cat meat case from the PSSE/PIW Oświęcim area (researched by the Jagiellonian University) is clustered with the 3-city subcluster. It is quite strange. History of traveling to the seaside or the Lublin region at the end of May or cross-contamination of material from Gdańsk (because it also passed through their lab) should be done. As I noticed earlier, this Lublin region mixes well and epidemiologically it looks as if all paths came from it/through it (even Lublin cases are clustered with cases from the Nysa and Oder river basins). ## 2.4 DBScan of epidemic inference based on spatial, temproral and genetic distances

We take 19 samples from WOAH reference lab for detailed investigation.

 cats_ref_norm=norm_scale
 cats_ref_norm$z=c(1:57)
piwet_norm=norm_scale[c(1:19),]
piwet_norm[1,]=cats_ref_norm[3,]
piwet_norm[2,]=cats_ref_norm[4,]
piwet_norm[3,]=cats_ref_norm[7,]
piwet_norm[4,]=cats_ref_norm[8,]
piwet_norm[5,]=cats_ref_norm[11,]
piwet_norm[6,]=cats_ref_norm[17,]
piwet_norm[7,]=cats_ref_norm[19,]
piwet_norm[8,]=cats_ref_norm[20,]
piwet_norm[9,]=cats_ref_norm[23,]
piwet_norm[10,]=cats_ref_norm[24,]
piwet_norm[11,]=cats_ref_norm[28,]
piwet_norm[12,]=cats_ref_norm[5,]
piwet_norm[13,]=cats_ref_norm[6,]
piwet_norm[14,]=cats_ref_norm[1,]
piwet_norm[15,]=cats_ref_norm[2,]
piwet_norm[16,]=cats_ref_norm[12,]
piwet_norm[17,]=cats_ref_norm[9,]
piwet_norm[18,]=cats_ref_norm[10,]

dist=matrix(data=NA, nrow=19, ncol=19)
dist_t=matrix(data=NA, nrow=19, ncol=19)

for (i in 1:19)
for (j in 1:19)
{
  {
  dist[i,j]=sqrt((piwet_norm$lon[i]-piwet_norm$lon[j])^2+ (piwet_norm$lat[i]-piwet_norm$lat[j])^2)
  dist_t[i,j]=abs(piwet_norm$week[i]-piwet_norm$week[j])
  } 
  }
#impact of each dimensions to
klaster=p_piwet/2+2*dist+ 2*dist_t
  ds <- dbscan(klaster, eps = 5, minPts = 2)
ds
## DBSCAN clustering for 19 objects.
## Parameters: eps = 5, minPts = 2
## Using euclidean distances and borderpoints = TRUE
## The clustering contains 2 cluster(s) and 7 noise points.
## 
## 0 1 2 
## 7 8 4 
## 
## Available fields: cluster, eps, minPts, dist, borderPoints
pheatmap(klaster, display_numbers = T,labels_row = l)

We see that 1 Poznan case (upper left corner) is from a separated introduction of virus. One big Blocking cluster of Western Poland can be seen (with a one linking case from a Lublin region). Weakly connected cluster of Lublin cases is complementing the hierarchical matrix.

plot(newmap, xlim = c(15, 25), ylim = c(45, 55), asp = 1)
points(koty_sorter$lon[piwet_norm$z] +rnorm(19)/10, koty_sorter$lat[piwet_norm$z] +rnorm(19)/10, col = ds$cluster+1, pch=20, cex = 1)
title("DBscan genetics/2+2*geography+ 2*time")

klaster=p_piwet/3+3*dist+ 1*dist_t
  ds <- dbscan(klaster, eps = 6, minPts = 2)
ds
## DBSCAN clustering for 19 objects.
## Parameters: eps = 6, minPts = 2
## Using euclidean distances and borderpoints = TRUE
## The clustering contains 2 cluster(s) and 1 noise points.
## 
##  0  1  2 
##  1 11  7 
## 
## Available fields: cluster, eps, minPts, dist, borderPoints
plot(newmap, xlim = c(15, 25), ylim = c(45, 55), asp = 1)
points(koty_sorter$lon[piwet_norm$z] +rnorm(19)/10, koty_sorter$lat[piwet_norm$z] +rnorm(19)/10, col = ds$cluster+1, pch=20, cex = 1)
title("DBscan genetics/3+3*geography+ 1*time")

There are 3 clear clusters (sensitive to choice of weights in inference model). First is a single case in Poznan, second western Poland and third Lublin area. Second and Third clusters are interlinked.

3. Conclusions of data triangulation

Some of results of different methods suggest opposite routes. One need to remember that change in epizootic dynamics may be driven by animal owners awareness and this kind of correction my cover biological processes. Due to the lack of official data accurate with time and place of material collection, but only the testing time since 22.06, it is worth looking at unofficial data before. Thanks to participatory study and media analytics, we could reconstruct the most probable version of the first phases of the outbreak.

3.1 Main outcomes

If one can take into account possible biased of each method then, data triangulation suggest that most likely (according to the Google and media analysis, reference lab results and participatory epidemiology project with use only of spatio-temporal information):

  • the first cases A/H5N1 were probably in cats in second half of May in regions bordering with Ukraine (possible meat and wild life sources). This information could be missing without retrospective media analysis and participatory reporting;

  • the biggest burden of disease was in Pomerania (in particular in Gdansk) in the middle of June;

  • positivity ratio is the highest in Western Poland, which suggest that the last outbreak clusters were in end June/beginning of July in Western Poland;

  • the positive cases are clearly forming chains (in contradiction to negative cases and to participatory reports) which support environmental source hypothesizes;

  • most of positive cases are both on bird migration paths and close to high abundance of nesting sites of water birds;

  • poultry meat contamination as additional factor of spread would explain some additional transmissions (i.e. cat case in Łodź or cannot be easily explained by environmental exposure).

3.2 Water bird migration and nesting

We suggest environmental or water birds active monitoring in Lubelskie province (probably there is no more viral pressure as the peak of cases was in May) and Gdańsk area (probably there is no more viral pressure as the cases were in June).

Some migration routes in middle/late Spring (before nesting season) are more less known. As well as geographical abundance of species or nesting site can be found in birds atlas or citizens participatory ornithological projects.

par(mfrow=c(1,3))
knitr::include_graphics("breeding_all.png", dpi=100)

knitr::include_graphics("migration.png", dpi=100)

knitr::include_graphics("breeding_duck.png", dpi=100)

knitr::include_graphics("risks.png", dpi=100)

3.2 Meat food chains

We propose an analysis of possible contamination of meat/eggs imported from Ukraine or from farms/slaughterhouses in Lubelskie province. Source of meat positive sample from Oswiecim region should be determinated by backward tracing.

Example of chains can be easily found in the Internet and are available in ministry of agriculture registry.

par(mfrow=c(1,2))
knitr::include_graphics("dino.jpg", dpi=100)

knitr::include_graphics("lidl.jpg", dpi=100)

3.3 Extended epidemiolpogical analysis

Data from meat samples should be included into consideration. Genetic distance between strains could help in clustering analysis and source detection.

Travel history of the caracal case from Łódź must be determined. <update 03.08> This case was wrongly geocoded (there in error in article https://www.eurosurveillance.org/content/10.2807/1560-7917.ES.2023.28.31.2300366) and it was in Ełk indeed.

Possible excess mortality rates from veterinary clinics in May on Polish Ukrainian border should be analyzed.

Regression analysis for risk factors among positive and negative results should be performed. The study will report on the spatio-temporal dynamics and the factors that might have influenced A/H5N1 in cats that will be estimated in Poland.

Genetic distance to AI strains found in Poland in poultry and wild birds in April/May should be calculated.

We speculate, that virus found in meant in Oswiecim county (meat bought in beginning of June: EPI_ISL_17959737) should be close to (if our hypothesis of introduction thought Lublin/border area is valid) early cases from Lublin (EPI_ISL_17951051, EPI_ISL_17951052, EPI_ISL_17951053). <update 03.08> The cases from Oswiecim (near KrakóW) are genetically closed to some Lublin and Gdańsk cases.