A novel method for inference of acyclic chemical compounds with bounded branch-height based on artificial neural networks and integer programming

Analysis of chemical graphs is becoming a major research topic in computational molecular biology due to its potential applications to drug design. One of the major approaches in such a study is inverse quantitative structure activity/property relationship (inverse QSAR/QSPR) analysis, which is to infer chemical structures from given chemical activities/properties. Recently, a novel two-phase framework has been proposed for inverse QSAR/QSPR, where in the first phase an artificial neural network (ANN) is used to construct a prediction function. In the second phase, a mixed integer linear program (MILP) formulated on the trained ANN and a graph search algorithm are used to infer desired chemical structures. The framework has been applied to the case of chemical compounds with cycle index up to 2 so far. The computational results conducted on instances with n non-hydrogen atoms show that a feature vector can be inferred by solving an MILP for up to \(n=40\), whereas graphs can be enumerated for up to \(n=15\). When applied to the case of chemical acyclic graphs, the maximum computable diameter of a chemical structure was up to 8. In this paper, we introduce a new characterization of graph structure, called “branch-height” based on which a new MILP formulation and a new graph search algorithm are designed for chemical acyclic graphs. The results of computational experiments using such chemical properties as octanol/water partition coefficient, boiling point and heat of combustion suggest that the proposed method can infer chemical acyclic graphs with around \(n=50\) and diameter 30.

Source code of the implementation of our algorithm is freely available from https://github.com/ku-dml/mol-infer.

ANN:

Artificial neural network

MILP:

Mixed integer linear programming

Not applicable.

This research was supported, in part, by Japan Society for the Promotion of Science, Japan, under Grant #18H04113.

Authors and Affiliations

1. Department of Applied Mathematics and Physics, Kyoto University, Yoshida Honmachi, Sakyo, Kyoto, 606-8501, Japan

   Naveed Ahmed Azam, Jianshen Zhu, Yanming Sun, Yu Shi, Aleksandar Shurbevski & Hiroshi Nagamochi

2. Graduate School of Advanced Integrated Studies in Human Survivability, Kyoto University, Yoshida Nakaadachi-cho, Sakyo, Kyoto, 606-8306, Japan

   Liang Zhao

3. Bioinformatics Center, Institute for Chemical Research, Kyoto University, Gokasho, Uji, 611-0011, Japan

   Tatsuya Akutsu

Contributions

Conceptualization, HN and TA; methodology, HN; software, NAA, JZ, YS, YS, AS and L.; validation, NAA, JZ, AS and HN; formal analysis, HN; data resources, AS, LZ, HN and TA; writing—original draft preparation, HN; writing—review and editing, NAA, AS and TA; project administration, HN; funding acquisition, TA. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Hiroshi Nagamochi or Tatsuya Akutsu.

Competing interests

The authors declare that they have no competing interests.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Appendix A: Statistical features of molecular structures

We observe the following features of the graph-theoretical structure of chemical graphs registered in the chemical database PubChem. Let \(\text {DB}^{(\le n)}\) denote the set of chemical graphs with at most n non-hydrogen atoms that are registered in chemical database PubChem (downloaded a copy on March 21, 2019). The cycle index (or rank) of a chemical graph \(G=(H=(V,E),\alpha ,\beta )\) is defined to be \(|E|-(|V|-1)\) (i.e., the minimum number of edges to be removed to make the graph H acyclic). We call a chemical graph a rank-r chemical graph if the rank of the graph is r. The core of a chemical cyclic graph G is defined to be the induced subgraph \(G'\) of G such that \(G'\) consists of vertices in a cycle or vertices in a path joining two cycles. A vertex in the core (not in the core) is called a core vertex (resp., a non-core vertex). The edges not in the core of a chemical cyclic graph G form a collection of trees T, which we call a non-core tree. Each non-core tree contains exactly one core vertex and is regarded as a tree rooted at the core vertex. The k-branch height of a chemical cyclic graph G is defined to be the maximum of k-branch heights over all non-core trees.

Let \(\rho _r\) (%) denote the ratio of the number of chemical graphs with rank at most \(r\in [0,4]\) to the number of all chemical graphs in PubChem. See Table 11.

Let \(\rho _0^{(d)}\) (%) denote the ratio of the number of chemical graphs in \(\text {DB}^{(\le 100)}\) such that the maximum degree is at most \(d\in [3,4]\) to the number of all chemical graphs in \(\text {DB}^{(\le 100)}\). Let \(\rho _r^{(d)}\) (%), \(r\in [1,4]\) denote the ratio of the number of rank-r chemical graphs in \(\text {DB}^{(\le 100)}\) such that the maximum degree of a non-core vertex is at most \(d\in [3,4]\) to the number of all rank-r chemical graphs in \(\text {DB}^{(\le 100)}\). See Table 12.

Let \(\rho _r(k,h)\) (%), \(r\in [0,4]\), \(k=2\), \(h\in [1,2]\) denote the ratio of the number of rank-r chemical graphs in \(\text {DB}^{(\le 50)}\) such that the k-branch height is at most h to the number of all rank-r chemical graphs in \(\text {DB}^{(\le 50)}\). See Table 13. We see that most chemical graphs G with at most 50 non-hydrogen atoms satisfy \(\text {bh}_2(G)\le 2\).

We show the distribution of 2-branch height over alkans C\(_n\)H\(_{2n+2}\). Let \(\text {Aln}(n)\) denote the set of all alkans with n carbon atoms, where \(|\text {Aln}(25)|=36,797,588\). Let \(\rho _{\text {Aln}}(2,h)\) (%), \(h\in [1,4]\) denote the ratio of the number of alkans in \(\text {Aln}(25)\) such that the 2-branch height is at most h to the number of alkans in \(\text {Aln}(25)\). See Table 14.

Let \(\rho _{\text {2bt}}(\delta )\) denote the ratio of the number of acyclic chemical graphs in \(\text {DB}^{(\le 50)}\) such that the degree of the root of the 2-branch-tree is \(\delta \in [1,4]\) to the number of all acyclic chemical graphs in \(\text {DB}^{(\le 50)}\). See Table 15.

Among the 2-fringe-trees T of all acyclic chemical graphs in \(\text {DB}^{(\le 100)}\), over \(90\%\) of them satisfy \(n\le 2d+2\) for the number \(n=|V(T)|\) of non-hydrogen atoms in a 2-fringe-tree T and the number d of non-hydrogen atoms adjacent to the root in T.

Let \({\mathcal {FT}}_{0,2}\) denote the set of all 2-fringe-trees that appear in an acyclic chemical graph in \(\text {DB}^{(\le 100)}\), and \({\mathcal {FT}}_{0,2}^{(\delta )}\), \(\delta \in [1,3]\) denote the set of all 2-fringe-trees \(T\in {\mathcal {FT}}_{0,2}\) that have \(\delta \) children (i.e., the degree of the root is \(\delta \)). Let \(\rho _{2\delta +2}^{(\delta )}\) (%) denote the ratio of the number of 2-fringe-trees in \({\mathcal {FT}}_{0,2}^{(\delta )}\) that have at most \(2\delta +2\) vertices to the number of 2-fringe-trees in \({\mathcal {FT}}_{0,2}^{(\delta )}\). See Table 16.

Appendix B: Formulating an MILP based on scheme graphs

This section shows how to formulate an MILP based on a scheme graph.

Scheme graphs

Let \(t^*\), \(s^*\), and \(c^*\), be integers such that

• \(t^*= n^* - (\text {bh}^* -1) - (k^*+1){\text {bl}}^*\);

• \(s^*=a(b^c-1)/(b-1)+1\) for \(a=d_\text {max}\), \(b=d_\text {max}-1\) and \(c=\text {bh}^*\); and

• \(c^*=s^*-1\).

Let a scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*)\) consist of a tree \(T_B\), a path \(P_{t^*}\), a set \(\{S_s\mid s\in [1,s^*]\}\) of trees, a set \(\{T_t\mid t\in [1,t^*]\}\) of trees, and a set of directed edges between \(T_B\) and \(P_{t^*}\) so that an acyclic graph \(H\in {\mathcal {H}}(n^*, d_\text {max}, \text {dia}^*, k^*, \text {bh}^*, {\text {bl}}^*)\) will be constructed in the following way:

1. (i)

   The \(k^*\)-branch-tree of H will be chosen as a subtree of \(T_B=(V_B,E_B)\);

2. (ii)

   Each \(k^*\)-fringe-tree rooted at a vertex \(u_s\in V(T_B)\) of H will be chosen as a subtree of \(S_s\);

3. (iii)

   Each \(k^*\)-branch-path of H (except for its end-vertices) will be chosen as a subpath of \(P_{t^*}\) or as an edge in \(T_B\);

4. (iv)

   Each \(k^*\)-fringe-tree rooted at a vertex \(v_t\in V(P_{t^*})\) of H will be chosen as a subtree of \(T_t\); and

5. (v)

   An edge (u, v) directed from \(T_B\) to \(P_{t^*}\) will be selected as an initial edge of a \(k^*\)-branch-path of H and an edge (v, u) directed from \(P_{t^*}\) to \(T_B\) will be selected as an ending edge of a \(k^*\)-branch-path of H.

More formally, each component of a scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*)\) is defined as follows.

1. (i)

   \(T_B =(V_B=\{u_1,u_2,\ldots ,u_{s^*}\},\) \(E_B=\{a_1,a_2,\ldots ,a_{c^*}\})\), called a base-tree is a tree rooted at a vertex \(u_1\) that is isomorphic to the rooted tree \(T(d_\text {max}, d_\text {max}- 1, \text {bh}^*)\). Regard \(T_B\) as an ordered tree by introducing a total order for each set of siblings and call the first (resp., last) child in a set of siblings the leftmost (resp. rightmost) child, which defines the leftmost (rightmost) path from the root \(u_1\) to a leaf in \(T_B\), as illustrated in Fig. 4a.

   For each vertex \(u_s\in V_B\), let \(E_B(s)\) denote the set of indices i of edges \(a(i)\in E_B\) incident to \(u_s\) and \(\text {Cld}_B(s)\) denote the set of indices i of children \(u_i\in V_B\) of \(u_s\) in the tree \(T_B\).

   For each integer \(d\in [0,k^*]\), let \(V_B(d)\) denote the set of indices s of vertices \(u_s\in V_B\) whose depth is d in the tree \(T_B\), where \(V_B(\text {bh}^*)\) is the set of indices s of leaves \(u_s\) of \(T_B\).

   Regard each edge \(a_i\in E_B\) as a directed edge \((u_s,u_{s'})\) from one end-vertex \(u_s\) of \(a_i\) to the other end-vertex \(u_{s'}\) of \(a_i\) such that \(s=\text {prt}(s')\) (i.e., \(u_s\) is the parent of \(u_{s'}\)), where \(\text {head}(i)\) and \(\text {tail}(i)\) denote the head \(u_{s'}\) and tail \(u_s\) of edge \(a_i\in E_B\), respectively.

   For each index \(s\in [1,s^*]\), let \(E_B^+(s)\) (resp., \(E_B^-(s)\)) denote the set of indices i of edges \(a_i\in E_B\) such that the tail (resp., head) of \(a_i\) is vertex \(u_{s}\).

   Let \(L_B\) denote the set of indices of leaves of \(T_B\), and \(s^\text {left}\) (resp., \(s^\text {right}\)) denote the index \(s\in L_B\) of the leaf \(u_s\) at which the leftmost (resp., rightmost) path from the root ends.

   For each leaf \(u_s\), \(s\in L_B\), let \(V_{B,s}\) (resp., \(E_{B,s}\)) denote the set of indices s of non-root vertices \(u_s\) (resp., indices i of edges \(a(i)\in E_B\)) along the path from the root to the leaf \(u_s\) in the tree \(T_B\).

   For the example of a base-tree \(T_B\) with \(\text {bh}^*=2\) in Fig. 4, it holds that \(L_B=\{5,6,7,8,9,10\}\), \(s^\text {left}=5\), \(s^\text {right}=10\), \(E_{B,s^\text {left}}=\{1,4\}\) and \(V_{B,s^\text {left}}=\{2,5\}\).

2. (ii)

   \(S_s\), \(s\in [1,s^*]\) is a tree rooted at vertex \(u_s\in V_B\) in \(T_B\) that is isomorphic to the rooted tree \(T(d_\text {max}-1, d_\text {max}-1, k^*)\), as illustrated in Fig. 4b. Let \(u_{s,i}\) and \(e'_{s,i}\) denote the vertex and edge in \(S_s\) that correspond to the i-th vertex and the i-th edge in \(T(d_\text {max}-1, d_\text {max}-1, k^*)\), respectively. Regard each edge \(e'_{s,i}\) as a directed edge \((u_{s,\text {prt}(i)},u_{s,i})\). For this, each vertex \(u_s\in V_B\) is also denoted by \(u_{s,1}\).

3. (iii)

   \(P_{t^*}=(V_P=\{v_1\), \(v_2\), \(\ldots \), \(v_{t^*}\}\), \(E_P=\{e_2\), \(e_3\), \(\ldots \), \(e_{t^*}\})\), called a link-path with size \(t^*\) is a directed path from vertex \(v_1\) to vertex \(v_{t^*}\), as illustrated in Fig. 4a. Each edge \(e_t\in E_P\) is directed from vertex \(v_{t-1}\) to vertex \(v_t\).

4. (iv)

   \(T_t\), \(t\in [1,t^*]\) is a tree rooted at vertex \(v_t\) in \(P_{t^*}\) that is isomorphic to the rooted tree \(T(d_\text {max}-2, d_\text {max}-1, k^*)\), as illustrated in Fig. 4c. Let \(v_{t,i}\) and \(e_{t,i}\) denote the vertex and edge in \(T_t\) that correspond to the i-th vertex and the i-th edge in \(T(d_\text {max}-2, d_\text {max}-1, k^*)\), respectively. Regard each edge \(e_{t,i}\) as a directed edge \((v_{t,\text {prt}(i)},u_{t,i})\). For this, each vertex \(v_t\in V_P\) is also denoted by \(v_{t,1}\).

5. (v)

   For every pair (s, t) with \(s\in [1,s^*]\) and \(t\in [1,t^*]\), join vertices \(u_{s}\) and \(v_{t}\) with directed edges \((u_{s},v_{t})\) and \((v_{t},u_{s})\), as illustrated in Fig. 4a.

We explain the basic idea of an MILP in Theorem 2. The MILP mainly consists of the following three types of constraints.

1. C1.

   Constraints for selecting an acyclic graph H as a subgraph of the scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*);\)

2. C2.

   Constraints for assigning chemical elements to vertices and multiplicity to edges to determine a chemical graph \(G=(H,\alpha ,\beta )\); and

3. C3.

   Constraints for computing descriptors from the selected acyclic chemical graph G.

In the constraints of C1, more formally we prepare the following.

1. (i)

   In the scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*)\), we prepare a binary variable u(s, 1) for each vertex \(u_s=u_{s,1}\in V_B\), \(s\in [1,s^*]\) so that vertex \(u_s=u_{s,1}\) becomes a \(k^*\)-branch of a selected graph H if and only if \(u(s,1)=1\). The subgraph of the base-tree \(T_B\) that consists of vertices \(u_s=u_{s,1}\) with \(u(s,1)=1\) will be the \(k^*\)-branch-tree of the graph H. We also prepare a binary variable a(i), \(i\in [1,c^*]\) for each edge \(a_i\in E_B\), where \(c^*=s^*-1\). For a pair of a vertex \(u_{s,1}\) and a child \(u_{s',1}\) of \(u_{s,1}\) such that \(u(s,1)=u(s',1)=1\), either the edge \(a_i=(u_{s,1},u_{s',1})\) is used in the selected graph H (when \(a(i)=1\)) or a path \(P_i=(u_{s,1},v_{t',1},v_{t'+1,1},\ldots ,v_{t'',1},u_{s',1})\) from vertex \(u_{s,1}\) to vertex \(u_{s',1}\) is constructed in H with an edge \((u_{s,1},v_{t',1})\), a subpath \((v_{t',1},v_{t'+1,1},\ldots ,v_{t'',1})\) of the link-path \(P_{t^*}\) and an edge \((v_{t'',1},u_{s',1})\) (when \(a(i)=0\)). For example, vertices \(u_{1,1}\) and \(u_{2,1}\) are connected by a path \(P_1=(u_{1,1},v_{1,1},v_{2,1},u_{2,1})\) in the selected graph \(H'\) in Fig. 5c.

2. (ii)

   Let

   • \(n_\text {tree}^{\text {S}} =1+(d_\text {max}-1)((d_\text {max}-1)^{k^*}-1)/(d_\text {max}-2),\)

   • \(n_\text {tree}^{\text {T}} =1+(d_\text {max}-2)((d_\text {max}-1)^{k^*}-1)/(d_\text {max}-2),\)

   where \(n_\text {tree}^{\text {S}} \) (resp., \(n_\text {tree}^{\text {T}}\)) is the number of vertices in the rooted tree \(T(d_\text {max}-1, d_\text {max}-1, k^*)\) (resp., \(T(d_\text {max}-2, d_\text {max}-1, k^*)\)). In each tree \(S_s\), \(s\in [1,s^*]\) (resp., \(T_t\), \(t\in [1,t^*]\)) in the scheme graph, we prepare a binary variable u(s, i) (resp., v(t, i)) for each vertex \(u_{s,i}\), \(i\in [2, n_\text {tree}^{\text {S}}]\) (resp., \(v_{t,i}\), \(i\in [2, n_\text {tree}^{\text {T}}]\)) so that \(u(s,i)=1\) (resp., \(v(t,i)=1\)) means that the corresponding vertex \(u_{s,i}\) (resp., \(v_{t,i}\)) is used as a vertex in a selected graph H. The (non-empty) subgraph of a tree \(S_s\) (resp., \(T_t\)) that consists of vertices \(u_{s,i}\) with \(u(s,i)=1\) (resp., \(v_{t,i}\) with \(v(t,i)=1\)) will be a \(k^*\)-fringe-tree of a selected graph H.

3. (iii)

   In the link-path \(P_{t^*}\), we prepare a binary variable e(t), \(t\in [2,t^*]\) for each edge \(e_{t,1}=(v_{t-1,1}, v_{t,1})\in E_P\) so that \(e(t)=1\) if and only if edge \(e_{t,1}\) is used in some path \(P_i=(u_{s,1},v_{t',1},v_{t'+1,1},\ldots ,v_{t'',1},u_{s',1})\) constructed in (i).

4. (iv)

   For each pair (s, t) of \(s\in [1,s^*]\) and \(t\in [1,t^*]\), we prepare a binary variable e(s, t) (resp., e(t, s)) so that \(e(s,t')=1\) (resp., \(e(t'',s)=1\)) if and only if directed edge \((u_{s,1},v_{t',1})\) (resp., \((v_{t'',1},u_{s,1})\)) is used as the first edge (resp., last edge) of some path \(P_i=(u_{s,1},v_{t',1},v_{t'+1,1},\ldots , v_{t'',1},u_{s',1})\) constructed in (i).

Based on these, we include constraints with some more additional variables so that a selected subgraph H is a connected acyclic graph. See constraints (12) to (32) in Appendix C for the details.

In the constraints of C2, we prepare an integer variable \({\widetilde{\alpha }}(u)\) for each vertex u in the scheme graph that represents the chemical element \(\alpha (u)\in \Lambda \) if u is in a selected graph H (or \({\widetilde{\alpha }}(u)=0\) otherwise) and an integer variable \({\widetilde{\beta }}(e)\in [0,3]\) (resp., \({\widehat{\beta }}(e)\in [0,3]\)) for each edge e (resp., \(e=e(s,t)\) or e(t, s), \(s\in [1,s^*]\), \(t\in [1,t^*]\)) in the scheme graph that represents the multiplicity \(\beta (e)\in [1,3]\) if e is in a selected graph H (or \({\widetilde{\beta }}(e)\) or \({\widehat{\beta }}(e)\) takes 0 otherwise). This determines a chemical graph \(G=(H,\alpha ,\beta )\). Also we include constraints for a selected chemical graph G to satisfy the valence condition \((\alpha (u),\alpha (v),\beta (uv))\in \Gamma \) for each edge \(uv\in E\). See constraints (33) to (47) in Appendix C for the details.

In the constraints of C3, we introduce a variable for each descriptor and constraints with some more variables to compute the value of each descriptor in f(G) for a selected chemical graph G. See constraints (48) to (75) in Appendix C for the details.

Appendix C: All constraints in an MILP formulation for chemical acyclic graphs

To formulate an MILP that represents a chemical graph, we distinguish a tuple \((\mathtt{a,b},m)\) from a tuple \((\mathtt{b,a},m)\). For a tuple \(\gamma =(\mathtt{a,b},m)\in \Lambda \times \Lambda \times \{1,2,3\}\), let \(\overline{\gamma }\) denote the tuple \((\mathtt{b, a}, m)\). Let \(\Gamma _{<}\triangleq \{\overline{\gamma }\mid \gamma \in \Gamma _{>}\}\). We call a tuple \(\gamma =(\mathtt{a, b}, m) \in \Lambda \times \Lambda \times \{1, 2, 3\}\) proper if \(m \le \min \{ {\text {val}}({\mathtt{a}}), {\text {val}}(\mathtt{b})\}\) and \(m \le \max \{ {\text {val}}({\mathtt{a}}), {\text {val}}(\mathtt{b}) \}-1\), where the latter is assumed because otherwise G must consist of two atoms of \(\mathtt{a=b}\). Assume that each tuple \(\gamma \in \Gamma \) is proper. Let \(\epsilon \) be a fictitious chemical element that represents null, call a tuple \((\mathtt{a,b},0)\) with \(\mathtt{a,b}\in \Lambda \cup \{\epsilon \}\) fictitious, and define \(\Gamma _0\) to be the set of all fictitious tuples; i.e., \(\Gamma _0=\{(\mathtt{a,b},0) \mid \mathtt{a,b}\in \Lambda \cup \{\epsilon \}\}\). To represent chemical elements \(\mathtt{e}\in \Lambda \cup \{\epsilon \}\cup \Gamma \) in an MILP, we encode these elements \(\mathtt{e}\) into some integers denoted by \([\mathtt{e}]\). Assume that, for each element \({\mathtt{a}}\in \Lambda \), \([{\mathtt{a}}]\) is a positive integer and that \([\epsilon ]=0\).

Upper and lower bounds on descriptors

In our formulation of an MILP for inferring a vector \(x^*\) in Stage 4, we fix the following descriptors as specified constants: the number n(G) of vertices, the diameter \(\text {dia}(G)\), and the number \({\text {bl}}_{k^*}(G)\) of leaf \(k^*\)-leaf branches, which are set to be given integers \(n^*\), \(\text {dia}^*\), and \({\text {bl}}^*\), respectively. For each of the other descriptors, we specify a lower bound \(\text {LB}\) and an upper bound \(\text {UB}\) on the value so that the descriptor takes a value from the range between \(\text {LB}\) and \(\text {UB}\).

constants

• \(n^*\ge 5\): the size n(G) of G;

• \(\text {LB}_\text {dg}^\text {t}(i), \text {UB}_\text {dg}^\text {t}(i)\in [0,n^*],\) \(i\in [1,4], \text {t}\in \{{\text {in}},\text {ex}\}\): lower and upper bounds on the number \(\text {dg}_i^\text {t}(G)\) of \(k^*\)-internal/\(k^*\)-external vertices of degree i in G;

• \(\text {LB}_\text {ce}^\text {t}({\mathtt{a}}), \text {UB}_\text {ce}^\text {t}({\mathtt{a}})\in [0,n^*]\), \({\mathtt{a}}\in \Lambda , \text {t}\in \{{\text {in}},\text {ex}\}\): lower and upper bounds on the number \(\text {ce}_{\mathtt{a}}^\text {t}(G)\) of \(k^*\)-internal/\(k^*\)-external vertices v with \(\alpha (v)={\mathtt{a}}\) in G;

• \(\text {LB}_\text {bd}^\text {t}(m)\), \(\text {UB}_\text {bd}^\text {t}(m)\in [0,n^*-1],\) \(m\in [2,3],\)\(\text {t}\in \{{\text {in}},\text {ex}\}\): lower and upper bounds on the number \(\text {bd}_m^\text {t}(G)\) of \(k^*\)-internal/\(k^*\)-external edges e with \(\beta (e)=m\) in G;

• \(\text {LB}_\text {ac}^\text {t}(\gamma ),\) \(\text {UB}_\text {ac}^\text {t}(\gamma )\in [0,n^*-1],\) \(\text {t}\in \{{\text {in}},\text {ex}\},\) \(\gamma \in \Gamma _{<}\cup \Gamma _{=}\): lower and upper bounds on the number \( \text {ac}_{\gamma }^\text {t}(G)\) of \(k^*\)-internal/\(k^*\)-external edges e with adjacency-configuration \(\gamma \) in G;

• \(\text {LB}_\text {bc}^\text {t}(\mu ), \text {UB}_\text {bc}^\text {t}(\mu )\in [0,n^*-1], \text {t}\in \{{\text {in}},\text {ex}\}\), \(\mu \in \text {Bc}\): lower and upper bounds on the number \(\text {bc}_{\mu }^\text {t}(G)\) of \(k^*\)-internal/\(k^*\)-external edges e with bond-configuration \(\mu \) in G;

variables x for descriptors

• \(\text {dg}^{\text {in}}(i), \text {dg}^\text {ex}(i)\in [0,n^*]\), \(i\in [1,4]\): \(\text {dg}^{\text {in}}(i)\) (resp., \(\text {dg}^\text {ex}(i)\)) represents \(\text {dg}_i^{\text {in}}(G)\) (resp., \(\text {dg}_i^\text {ex}(G)\));

• \(\text {ce}^{\text {in}}({\mathtt{a}}), \text {ce}^\text {ex}({\mathtt{a}})\in [0,n^*]\), \({\mathtt{a}}\in \Lambda \): \(\text {ce}^{\text {in}}({\mathtt{a}})\) (resp., \(\text {ce}^\text {ex}({\mathtt{a}})\)) represents \(\text {ce}_{\mathtt{a}}^{\text {in}}(G)\) (resp., \(\text {ce}_{\mathtt{a}}^\text {ex}(G)\));

• \(\text {bd}^{\text {in}}(m), \text {bd}^\text {ex}(m)\in [0,2n^*]\), \(m\in [1,3]\): \(\text {bd}^{\text {in}}(m)\) (resp., \(\text {bd}^\text {ex}(m)\)) represents \(\text {bd}_m^{\text {in}}(G)\) (resp., \(\text {bd}_m^\text {ex}(G)\));

• \(\text {ac}^{\text {in}}(\gamma ), \text {ac}^\text {ex}(\gamma )\in [0,n^*]\), \(\gamma \in \Gamma _{<}\cup \Gamma _{=}\): \(\text {ac}^{\text {in}}(\gamma )\) (resp., \(\text {ac}^\text {ex}(\gamma )\)) represents \( \text {ac}_{\gamma }^{\text {in}}(G)\) (resp., \( \text {ac}_{\gamma }^\text {ex}(G)\));

• \(\text {bc}^{\text {in}}(\mu ), \text {bc}^\text {ex}(\mu )\in [0,n^*-1]\), \(\mu \in \text {Bc}\): \(\text {bc}^{\text {in}}(\mu )\) (resp., \(\text {bc}^\text {ex}(\mu )\)) represents \(\text {bc}_{\mu }^{\text {in}}(G)\) (resp., \(\text {bc}_{\mu }^\text {ex}(G)\));

constraints

We use the range-based method to define an applicability domain for our method. For this, we find the range (the minimum and maximum) of each descriptor over all relevant chemical compounds and represent each range as a set of linear constraints in the constraint set \({\mathcal {C}}_1\) of our MILP formulation. Recall that \(D_{\pi }\) stands for a set of chemical graphs used for constructing a prediction function. However, the number of examples in \(D_{\pi }\) may not be large enough to capture a general feature on the structure of chemical graphs. For this, we also use some data set from the whole set \(\text {DB}\) of chemical graphs in a database. Let \(\text {DB}_{{\mathcal {G}}}^{(i)}\) denote the set of chemical graphs \(G\in \text {DB}\cap {\mathcal {G}}\) such that \(n(G)=i\) for each integer \(i\ge 1\). Based on this, we assume that the given lower and upper bounds on the above descriptors satisfy the following. For each \(\text {t}\in \{{\text {in}},\text {ex}\}\),

Construction of scheme graph

We infer a subgraph H such that the maximum degree is \(d_\text {max}\in \{3,4\}\), \(n(H)=n^*\), \(\text {bh}_{k^*}(H)=\text {bh}^*\), and \({\text {bl}}_{k^*}(H)={\text {bl}}^*\). For this, we first construct the scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*)\). We then prepare a binary variable u(s, i) (resp., v(t, i)) for each vertex \(u_{s,i}\) in tree \(S_s\) (resp., \(v_{t,i}\) in tree \(T_t\)).

Recall that when the two end-vertices of edge \(a_i=(u_{s,1},u_{s',1}) \in E_B=\{a_1,a_2,\ldots ,a_{c^*}\}\) is connected in a selected subgraph H, either edge \(a_i\) is directly used in H or a path \(P_i=(u_{s,1},v_{t',1},v_{t'+1,1},\ldots ,\) \(v_{t'',1},u_{s',1})\) from \(u_{s,1}\) to \(u_{s',1}\) visiting some vertices in \(P_{t^*}\) is constructed in H. We regard the index i of each edge \(a_i\in E_B=\{a_1,a_2,\ldots ,a_{c^*}\}\) as the “color” of the edge, and define the color set of \(E_B\) to be \([1,c^*]\). To introduce necessary linear constraints that can construct such a path \(P_i\) properly in our MILP, we assign the color i to the vertices \(v_{t',1},v_{t'+1,1},\ldots ,\) \(v_{t'',1}\) in \(P_{t^*}\) when a path \(P_i=(u_{s,1},v_{t',1},v_{t'+1,1},\ldots ,\) \(v_{t'',1},u_{s',1})\) is used in H.

constants

Integers \(d_\text {max}\in \{3,4\}\), \(n^*\ge 3\), \(\text {dia}^*\ge 3\), \(k^*\ge 1\), \(\text {bh}^*\ge 1\) and \({\text {bl}}^*\ge 2\);

variables

• \(a(i)\in \{0,1\}\), \(i\in E_B\): a(i) represents edge \(a_i\in E_B\) (\(a(i)=1\), \(i\in E_B\)) (\(a(i)=1\) \(\Leftrightarrow \) edge \(a_i\) is used in H);

• \(e(s,t),e(t,s)\in \{0,1\}\), \(s\in [1,s^*]\), \(t\in [1,t^*]\): e(s, t) (resp., e(t, s)) represents direction \((u_{s,1}, v_{t,1})\) (resp., \((v_{t,1}, u_{s,1})\)), where \(e(s,t)=1\) (resp., \(e(t,s)=1\)) \(\Leftrightarrow \) edge \(u_{s,1},v_{t,1}\) is used in H and direction \((u_{s,1}, v_{t,1})\) (resp., \((v_{t,1}, u_{s,1})\)) is assigned to edge \(u_{s,1} v_{t,1}\);

• \(\chi (t)\in [0,c^*]\), \(t\in [1,t^*]\): \(\chi (t)\) represents the color \(c\in [0,c^*]\) assigned to vertex \(v_{t,1}\) (\(\chi (t)=c\) \(\Leftrightarrow \) vertex \(v_{t,1}\) is assigned color c, where \(\chi (t)=c=0\) iff \(v_{t,1}\) is not in H);

• \(\delta _{\text {clr}}(t,c)\in \{0,1\}\), \(t\in [1,t^*]\), \(c\in [0,c^*]\) (\(\delta _{\text {clr}}(t,c)=1\) \(\Leftrightarrow \) \(\chi (t)=c\));

• \(\text {clr}(c)\in [0,t^*]\), \(c\in [0,c^*]\): the number of vertices \(v_{t,i}\) with color c;

• \(\deg ^\text {b+}(s)\in [0,4]\), \(s\in [1,s^*]\): the out-degree of vertex \(u_{s,1}\) in the \(k^*\)-branch-subtree of H;

• \(\deg ^\text {b-}(s)\in [0,4]\), \(s\in [1,s^*]\): the in-degree of vertex \(u_{s,1}\) in the \(k^*\)-branch-subtree of H;

constraints

Selecting a subgraph

From the scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*)\), we select a subgraph H such that \(n(H)=n^*\), \(\text {dia}(H)=\text {dia}^*\), \(\text {bh}_{k^*}(H)=\text {bh}^*\), and \({\text {bl}}_{k^*}(H)={\text {bl}}^*\).

constants

• Integers \(d_\text {max}\in \{3,4\}\), \(n^*\ge 3\), \(\text {dia}^*\ge 3\), \(k^*\ge 1\), \(\text {bh}^*\ge 1\) and \({\text {bl}}^*\ge 2\);

• For each tree \(S_s=T(d_\text {max}-1, d_\text {max}-1, k^*)\), prepare

  • the set \(\text {Cld}_{\text {S}}(i)\) of the indices of children of a vertex \(v_i\);

  • the index \(\text {prt}(i)\) of the parent of a non-root vertex \(v_i\);

  • the set \(\text {Dsn}_S(d)\) of indices i of a vertex \(v_i\) whose depth is d;

  • a proper set \(P_{\text {prc}}(d_\text {max}-1, d_\text {max}-1, k^*)\) of index pairs,

  • where we denote \(P_{\text {prc}}(d_\text {max}-1, d_\text {max}-1, k^*)\) by \(P_{S,\text {prc}}\);

• For each tree \(T_t=T(d_\text {max}-2, d_\text {max}-1, k^*)\), prepare

  • the set \(\text {Cld}_{\text {T}}(i)\) of the indices of children of a vertex \(v_i\);

  • the index \(\text {prt}(i)\) of the parent of a non-root vertex \(v_i\);

  • a proper set \(P_{\text {prc}}(d_\text {max}-2, d_\text {max}-1, k^*)\) of index pairs,

  • where we denote \(P_{\text {prc}}(d_\text {max}-2, d_\text {max}-1, k^*)\) by \(P_{T,\text {prc}}\);

variables

• \(\sigma (s)\in \{0,1\}\), \(s\in [1,s^* ]\): (\(\sigma (s)=1\) \(\Leftrightarrow \) vertex \(u_{s,1}\) is a non-leaf \(k^*\)-branch or a root);

• \(u(s,i)\in \{0,1\}\), \(s\in [1,s^* ]\), \(i\in [1,n_\text {tree}^{\text {S}}]\): u(s, i) represents vertex \(u_{s,i}\) (\(u(s,i)=1\) \(\Leftrightarrow \) vertex \(u_{s,i}\) is used in H and edge \(e'_{s,i}\) \((i\ge 2)\) is used in H), (\(u(s,1)=1\) and \(\sigma (s)=0\) \(\Leftrightarrow \) vertex \(u_{s,1}\) is a leaf \(k^*\)-branch);

• \(v(t,i)\in \{0,1\}\), \(t\in [1,t^* ]\), \(i\in [1,n_\text {tree}^{\text {T}}]\): v(t, i) represents vertex \(v_{t,i}\) (\(v(t,i)=1\) \(\Leftrightarrow \) vertex \(v_{t,i}\) is used in H and edge \(e_{t,i}\) \((i\ge 2)\) is used in H);

• \(e(t)\in \{0,1\}\), \(t\in [1,t^*+1]\): e(t) represents edge \(e_{t,1}=v_{t-1,1} v_{t,1}\), where \(e_{1,1}\) and \(e_{t^*+1,1}\) are fictitious edges (\(e(t)=1\) \(\Leftrightarrow \) edge \(e_{t,1}\) is used in H);

constraints

Constraints (21) and (22) represent an extension of constraint (1) on the size of 2-fringe-trees to the case of a general branch-parameter \(k^*\).

Assigning multiplicity

We prepare an integer variable \({\widetilde{\beta }}(e)\) or \({\widehat{\beta }}(e)\) for each edge e in the scheme graph \(\text {SG}( d_\text {max}, k^*, \text {bh}^*, t^*)\) to denote the multiplicity of e in a selected graph H and include necessary constraints for the variables to satisfy in H.

constants

• Prepare functions \(\text {tail}\) and \(\text {head}\) such that \(a_{i}=(u_{\text {tail}(i)}, u_{\text {head}(i)})\in E_B\);

• Assume that each edge in a tree \(S_s\), \(s\in [1,s^*]\) (resp., \(T_t\), \(t\in [1,t^*]\)) is denoted by \(e'_{s,i}\) (resp., \(e_{t,i}\)) with the integer \(i\in [2,n_\text {tree}^{\text {S}}]\) of the head \(u_{s,i}\) (resp., \(v_{t,i}\)) of the edge;

variables

• \({\widetilde{\beta }}(i)\in [0,3]\), \(i\in [1,c^*]\): \({\widetilde{\beta }}(i)\) represents the multiplicity of edge \(a_{i}\), where \({\widetilde{\beta }}(i)=0\) if edge \(a_{i}\) is not in an inferred chemical graph G;

• \({\widetilde{\beta }}(p,i)\in [0,3]\), \(p\in [1,s^*+ t^*]\), \(i\in [2,n_\text {tree}^{\text {S}}]\): \({\widetilde{\beta }}(p,i)\) with \(p\le s^*\) (resp., \(p>s^*\)) represents the multiplicity of edge \(e'_{p,i}\) (resp., \(e_{p-s^*,i}\));

• \({\widetilde{\beta }}(t,1)\in [0,3]\), \(t\in [1, t^*+1]\): \({\widetilde{\beta }}(t,1)\) represents the multiplicity of edge \(e_{t,1}\);

• \({\widehat{\beta }}(s,t)\in [0,3]\), \(s\in [1,s^*]\), \(t\in [1,t^*]\): \({\widehat{\beta }}(s,t)\) represents the multiplicity of edge \(u_{s,1}v_{t,1}\);

constraints

Assigning chemical elements and valence condition

We include constraints so that each vertex v in a selected graph H satisfies the valence condition; i.e., \(\beta (v)\le {\text {val}}(\alpha (v))\). With these constraints, a chemical acyclic graph \(G=(H,\alpha ,\beta )\) on a selected subgraph H will be constructed.

constants

• A set \(\Lambda \cup \{\epsilon \}\) of chemical elements, where \(\epsilon \) denotes null;

• A coding \([{\mathtt{a}}]\), \({\mathtt{a}}\in \Lambda \cup \{\epsilon \}\) such that \([\epsilon ]=0\); \([{\mathtt{a}}]\ge 1\), \({\mathtt{a}}\in \Lambda \); and \([{\mathtt{a}}]\ne [\mathtt{b}]\) if \({\mathtt{a}}\ne \mathtt{b}\); Let \([\Lambda ]\) and \([\Lambda \cup \{\epsilon \}]\) denote \(\{[{\mathtt{a}}]\mid {\mathtt{a}}\in \Lambda \}\) and \(\{[{\mathtt{a}}]\mid {\mathtt{a}}\in \Lambda \cup \{\epsilon \}\}\), respectively;

• A valence function: \({\text {val}}: \Lambda \rightarrow [1,4]\);

• Let \(E_B(s)\) denote the set of indices i of all edges \(a_i\in E_B\) adjacent to vertex \(u_{s,1}\) in \(T_B\).

variables

• \({\widetilde{\alpha }}(p,i)\in [\Lambda \cup \{\epsilon \}]\), \(p\in [1,s^*+ t^*]\), \(i\in [1,n_\text {tree}^{\text {S}}]\): \({\widetilde{\alpha }}(p,i)\) with \(p\le s^*\) (resp., \(p>s^*\)) represents \(\alpha (u_{p,i})\) (resp., \(\alpha (v_{p-s^*,i})\));

• \(\delta _{\alpha }(p,i,{\mathtt{a}})\in \{0,1\}\), \(p\in [1,s^*+ t^*]\), \(i\in [1,n_\text {tree}^{\text {S}}]\), \({\mathtt{a}}\in \Lambda \cup \{\epsilon \}\): \(\delta _{\alpha }(p,i,{\mathtt{a}})=1\) \(\Leftrightarrow \) \(\alpha (u_{p,i})={\mathtt{a}}\) for \(p\le s^*\) and \(\alpha (v_{p-s^*,i})={\mathtt{a}}\) for \(p> s^*\);

• \(\delta _{{\widetilde{\beta }}}(i,m)\in \{0,1\}\), \(p\in [1,s^*+ t^*]\), \(i\in [1,c^*]\), \(m\in [0,3]\): \(\delta _{{\widetilde{\beta }}}(i,m)=1\) \(\Leftrightarrow \) the multiplicity of edge \(a_i\) in an inferred chemical graph G is m;

• \(\delta _{{\widetilde{\beta }}}(p,i,m)\in \{0,1\}\), \(p\in [1,s^*+ t^*]\), \(i\in [2,n_\text {tree}^{\text {S}}]\), \(m\in [0,3]\): \(\delta _{{\widetilde{\beta }}}(p,i,m)=1\) \(\Leftrightarrow \) the multiplicity of edge \(e'_{p,i}\), \(p\le s^*\) (or \(e_{p-s^*,i}\), \(p>s^*\)) in G is m;

• \(\delta _{{\widetilde{\beta }}}(t,1,m)\in \{0,1\}\), \(t\in [1, t^*+1]\), \(m\in [0,3]\): \(\delta _{{\widetilde{\beta }}}(t,1,m)=1\) \(\Leftrightarrow \) the multiplicity of edge \(e_{t}\) in G is q;

• \(\delta _{{\widehat{\beta }}}(s,t,m)\in \{0,1\}\), \(s\in [1,s^*]\), \(t\in [1,t^*]\), \(m\in [0,3]\): \(\delta _{{\widehat{\beta }}}(s,t,m)=1\) \(\Leftrightarrow \) the multiplicity of edge \(u_{s,1}v_{t,1}\) in G is m;

constraints

Descriptors on mass, the numbers of elements and bonds

We include constraints to compute descriptors \(\overline{\text {ms}}(G)\), \(\text {ce}_{\mathtt{a}}(G)\) (\({\mathtt{a}}\in \Lambda )\), \(\text {bd}_m(G)\) (\(m\in [2,3]\)) and \(n_\mathtt{H}(G)\) according to the definitions in "Modeling of chemical compounds" section.

constants

• A function \(\text {mass}^*:\Lambda \rightarrow {\mathbb {Z}}\) (we let \(\text {mass}({\mathtt{a}})\) denote the observed mass of a chemical element \({\mathtt{a}}\in \Lambda \), and define \(\text {mass}^*({\mathtt{a}})=\lfloor 10\cdot \text {mass}({\mathtt{a}})\rfloor \));

variables

• \(\text {Mass}\in {\mathbb {Z}}\): \(\text {Mass}\) represents \(\sum _{v\in V} \text {mass}^*(\alpha (v))\);

• \(\text {bd}(m)\in [0,2n^*]\), \(m\in [1,3]\);

• \(\text {n}_\mathtt{H}\in [0,4n^*]\): the number \(n_\mathtt{H}(G)\) of hydrogen atoms to be included to G;

constraints

Descriptor for the Number of Specified Degree

We include constraints to compute descriptors \(\text {dg}_i(G)\) (\(i\in [1,4]\)) according to the definitions in "Modeling of chemical compounds" section. We also add constraints so that the maximum degree of a vertex in H is at most 3 (resp., equal to 4) when \(d_\text {max}=3\) (resp., \(d_\text {max}=4)\).

variables

• \(\deg (p,i)\in [0,4]\), \(p\in [1,s^*+ t^*]\), \(i\in [1,n_\text {tree}^{\text {S}}]\): \(\deg (p,i)\) represents \(\deg _H(u_{p,i})\) for \(p\le s^*\) or \(\deg _H(v_{p-s^*,i})\) for \(p> s^*\);

• \(\delta _{\deg }(p,i,d)\in \{0,1\}\), \(p\in [1,s^*+ t^*]\), \(i\in [1,n_\text {tree}^{\text {S}}]\), \(d\in [0,4]\): \(\delta _{\deg }(p,i,d)=1\) \(\Leftrightarrow \) \(\deg (p,i)=d\);

constraints

Descriptor for the number of adjacency-configurations

We include constraints to compute descriptors \(\text {ac}_{\gamma }(G)\) (\(\gamma =(\mathtt{a,b},m)\in \Gamma \)) according to the definitions in "Modeling of chemical compounds" section.

constants

• A set \(\Gamma =\Gamma _{<}\cup \Gamma _{=}\cup \Gamma _{>}\) of proper tuples \((\mathtt{a,b},m)\in \Lambda \times \Lambda \times [1,3]\);

• The set \(\Gamma _0=\{(\mathtt{a,b},0)\mid \mathtt{a,b}\in \Lambda \cup \{\epsilon \}\}\);

variables

• \(\delta _{\tau }(i,\gamma )\in \{0,1\}\), \(i\in [1,c^*]\), \(\gamma \in \Gamma \cup \Gamma _0 \): \(\delta _{\tau }(i, \gamma )=1\) \(\Leftrightarrow \) edge \(a_{i}\) is assigned tuple \(\gamma \); i.e., \(\gamma =({\widetilde{\alpha }}(\text {tail}(i),1), {\widetilde{\alpha }}(\text {head}(i),1), {\widetilde{\beta }}(i))\);

• \(\delta _{\tau }(t,1,\gamma )\in \{0,1\}\), \(t\in [2,t^*]\), \(\gamma \in \Gamma \cup \Gamma _0 \): \(\delta _{\tau }(t,1, \gamma )=1\) \(\Leftrightarrow \) edge \(e_{t,1}\) is assigned tuple \(\gamma \); i.e., \(\gamma =({\widetilde{\alpha }}(s^*+ t-1,1), {\widetilde{\alpha }}(s^*+ t,1), {\widetilde{\beta }}(t,1))\);

• \(\delta _{\tau }(p,i,\gamma )\in \{0,1\}\), \(p\in [1,s^*+ t^*]\), \(i\in [2,n_\text {tree}^{\text {S}}]\), \(\gamma \in \Gamma \cup \Gamma _0 \): \(\delta _{\tau }(p,i, \gamma )=1\) \(\Leftrightarrow \) edge \(e'_{p,i}\), \(p\le s^*\) (or \(e_{p-s^*,i}\), \(p> s^*\)) is assigned tuple \(\gamma \); i.e., \(\gamma =({\widetilde{\alpha }}(p,\text {prt}(i)), {\widetilde{\alpha }}(p,i), {\widetilde{\beta }}(p,i))\);

• \(\delta _{{\widehat{\tau }}}(s,t,\gamma )\in \{0,1\}\), \(s\in [1,s^*]\), \(t\in [1,t^*]\), \(\gamma \in \Gamma \cup \Gamma _0 \): \(\delta _{{\widehat{\tau }}}(s,t, \gamma )=1\) \(\Leftrightarrow \) edge \(u_{s,1}v_{t,1}\) is assigned tuple \(\gamma \); i.e., \(\gamma =({\widetilde{\alpha }}(s,1), {\widetilde{\alpha }}(s^*+ t,1), {\widehat{\beta }}(s,t))\);

constraints

Descriptor for bond-configuration

We include constraints to compute the descriptors for bond-configuration \(\text {bd}_{\mu }(G)\), \(\mu \in \text {Bc}\), according to the definition.

variables

• \(\text {bc}(\mu )\in [0,n^*-1]\), \(\mu \in \text {Bc}\);

• \(\delta _{\text {dc}}(i,d,d',m)\in \{0,1\}\), \(i\in [1,c^*]\), \(d,d'\in [0, 4]\), \(m\in [0,3]\): \(\delta _{\text {dc}}(i,d,d',m)=1\) \(\Leftrightarrow \) \(\deg _H(u_{{\text {tail}}(i)})=d\), \(\deg _H(u_{{\text {head}}(i)})=d'\) and \(\beta (a_i)=m\in [1,3]\) in G;

• \(\delta _{\text {dc}}(t,1,d,d',m)\in \{0,1\}\), \(t\in [2,t^*]\), \(d,d'\in [0, 4]\), \(m\in [0,3]\): \(\delta _{\text {dc}}(t,1,d,d',m)=1\) \(\Leftrightarrow \) \(\deg _H(v_{t-1,1})=d\), \(\deg _H(v_{t,1})=d'\) and \(\beta (e_{t,1})=m\in [1,3]\) in G;

• \(\delta _{\text {dc}}(p,i,d,d',m)\in \{0,1\}\), \(p\in [1,s^*+ t^*]\), \(i\in [2,n_\text {tree}^{\text {S}}]\), \(d,d'\in [0, 4]\), \(m\in [0,3]\): \(\delta _{\text {dc}}(p,i,d,d',m)=1\) \(\Leftrightarrow \) \(\deg _H(u_{p,\text {prt}(i)})=d\), \(\deg _H(u_{p,i})=d'\) and \(\beta (e'_{p,i})=m\in [1,3]\) for \(p\le s^*\) (or \(\deg _H(v_{p-s^*,\text {prt}(i)})=d\), \(\deg _H(v_{p-s^*,i})=d'\) and \(\beta (e_{p-s^*,i})=m\in [1,3]\) for \(p> s^*\)) in G;

• \(\delta _{\widehat{\text {dc}}}(s,t,d,d',m)\in \{0,1\}\), \(s\in [1,s^*]\), \(t\in [1,t^*]\), \(d,d'\in [0, 4]\), \(m\in [0,3]\): \(\delta _{\widehat{\text {dc}}}(s,t,d,d',1)=1\) \(\Leftrightarrow \) \(\deg _H(u_{s,1})=d\), \(\deg _H(v_{t,1})=d'\) and \(\beta (u_{s,1}v_{t,1})=m\in [1,3]\) in G;

constraints

Appendix D: Descriptions of new graph search algorithms

Multi-rooted trees and frequency vectors

For a finite set A of elements, let \({\mathbb {Z}}_+^A\) denote the set of functions \(\pmb {w}:A\rightarrow {\mathbb {Z}}_+\). A function \(\pmb {w}\in {\mathbb {Z}}_+^A\) is called a non-negative integer vector (or a vector) on A and the value \(\pmb {x}(a)\) for an element \(a\in A\) is called the entry of \(\pmb {x}\) for \(a\in A\). For a vector \(\pmb {w}\in {\mathbb {Z}}_+^A\) and an element \(a\in A\), let \(\pmb {w}+\pmb {1}_{a}\) (resp., \(\pmb {w}-\pmb {1}_{a}\)) denote the vector \(\pmb {w}'\) such that \(\pmb {w}'(a)=\pmb {w}(a)+1\) (resp., \(\pmb {w}'(a)=\pmb {w}(a)-1\)) and \(\pmb {w}'(b)=\pmb {w}(b)\) for the other elements \(b\in A\setminus \{a\}\). For a vector \(\pmb {w}\in {\mathbb {Z}}_+^A\) and a subset \(B\subseteq A\), let \(\pmb {w}_{[B]}\) denote the projection of \(\pmb {w}\) to B; i.e., \(\pmb {w}_{[B]}\in {\mathbb {Z}}_+^B\) such that \(\pmb {w}_{[B]}(b)=\pmb {w}(b)\), \(b\in B\).

Let \(\text {Bc}\) denote the set of tuples \(\mu =(d_1,d_2,k)\in [1,4]\times [1,4]\times [1,3]\) (bond-configuration) such that \(\max \{d_1,d_2\}+k\le 4\). For two tuples \(\mu =(d_1,d_2,k), \mu '=(d'_1,d'_2,k')\in \text {Bc}\), we write \(\mu \ge \mu '\) if

• \(\max \{d_1,d_2\}\ge \max \{d'_1,d'_2\}\), \(\min \{d_1,d_2\}\ge \min \{d'_1,d'_2\}\) and \(k\ge k'\),

and write \(\mu > \mu '\) if

• \(\mu \ge \mu '\) and \(\mu \ne \mu '\).

Let \(\text {Dg}=\{\text {dg}1, \text {dg}2, \text {dg}3, \text {dg}4\}\), where \(\text {dg}i\) denotes the number of vertices with degree i.

Henceforth we deal with vectors \(\pmb {w}\) that have their \(\pmb {w}_{\text {in}}\) and \(\pmb {w}_{\text {ex}}\) components, both \(\pmb {w}_{{\text {in}}}, \pmb {w}_{\text {ex}}\in {\mathbb {Z}}_+^{ {\Lambda \cup \Gamma \cup {\text {Bc}}\cup {\text {Dg}}}}\), and for convenience we write \(\pmb {w}= (\pmb {w}_{\text {in}}, \pmb {w}_\text {ex})\) in the sense of concatenation.

For a vector \(\pmb {x}= (\pmb {x}_{\text {in}}, \pmb {x}_\text {ex})\) with \(\pmb {x}_{\text {in}}, \pmb {x}_\text {ex}\in {\mathbb {Z}}_+^{ {\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}}\), let \({\mathcal {G}}(\pmb {x})\) denote the set of chemical acyclic graphs G whose 2-internal (resp., 2-external) vertices/edges are determined by the vector \(\pmb {x}_{\text {in}}\) (resp., \(\pmb {x}_\text {ex}\)); i.e., G satisfies the following:

• \(\text {ce}_{\mathtt{a}}^{\text {in}}(G)=\pmb {x}_{\text {in}}({\mathtt{a}})\) and \(\text {ce}_{\mathtt{a}}^\text {ex}(G)=\pmb {x}_\text {ex}({\mathtt{a}})\) for each chemical element \({\mathtt{a}}\in \Lambda \),

• \(\text {ac}_{\gamma }^{\text {in}}(G)= \pmb {x}_{\text {in}}(\gamma )\) and \(\text {ac}_{\gamma }^\text {ex}(G)=\pmb {x}_\text {ex}(\gamma )\) for each adjacency-configuration \(\gamma \in \Gamma \),

• \(\text {bc}_{\mu }^{\text {in}}(G)=\pmb {x}_{\text {in}}(\mu )\) and \(\text {bc}_{\mu }^\text {ex}(G)=\pmb {x}_\text {ex}(\mu )\) for each bond-configuration \(\mu \in \text {Bc}\),

• \(\text {dg}_i^{\text {in}}(G)=\pmb {x}_{\text {in}}(\text {dg}i)\) and \(\text {dg}_i^\text {ex}(G)=\pmb {x}_\text {ex}(\text {dg}i)\) for each degree \(\text {dg}i\in \text {Dg}\).

Throughout the section, let \(k^*=2\) be a branch-parameter, \(\pmb {x}^*= (\pmb {x}_{\text {in}}^*, \pmb {x}_\text {ex}^*)\) be a given feature vector with \(\pmb {x}_{\text {in}}^*, \pmb {x}_\text {ex}^*\in {\mathbb {Z}}_+^ {\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}\), and \(\text {dia}^*\) be an integer. We infer a chemical acyclic graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) such that \({\text {bl}}_2(G)\in [2,3]\) and the diameter of G is \(\text {dia}^*\), where \(n^*=\sum _{{\mathtt{a}}\in \Lambda }(\pmb {x}^*_{\text {in}}({\mathtt{a}})+\pmb {x}^*_\text {ex}({\mathtt{a}}) )\). Note that any other descriptors of \(G\in {\mathcal {G}}(\pmb {x}^*)\) can be determined by the entries of vector \(\pmb {x}^*\).

To infer a chemical acyclic graph \(G\in {\mathcal {G}}(\pmb {x}^*)\), we consider a connected subgraph T of G that consists of

Our method first generates a set \(\text {FT}\) of all possible rooted trees T that can be a 2-fringe-tree of a chemical graph \(G\in {\mathcal {G}}(\pmb {x}^*)\), and then extends the trees T by repeatedly appending a tree in \(\text {FT}\) until a chemical graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) is formed. In the extension, we actually manipulate the “frequency vectors” of trees defined below.

To specify which part of a given tree T plays the role of 2-internal vertices/edges or 2-external vertices/edges in a chemical graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) to be inferred, we designate at most three vertices \(r_1(T)\), \(r_2(T)\), and \(r_3(T)\), in T as terminals, and call T rooted (resp., bi-rooted and tri-rooted) if the number of terminals is one (resp., two and three). For a rooted tree (resp., bi- or tri-rooted tree) T, let \({\widetilde{V}}_{{\text {in}}}\) denote the set of vertices contained in a path between two terminals of T, \({\widetilde{E}}_{{\text {in}}}\) denote the set of edges in T between two vertices in \({\widetilde{V}}_{{\text {in}}}\), and define \({\widetilde{V}}_{\text {ex}}\triangleq V(T)\setminus {\widetilde{V}}_{{\text {in}}}\) and \({\widetilde{E}}_{\text {ex}}\triangleq E(T)\setminus {\widetilde{E}}_{{\text {in}}}\). For a bi- or tri-rooted tree T, define the backbone path \(P_T\) of T to be the path of T between vertices \(r_1(T)\) and \(r_2(T)\).

Given a chemical acyclic graph T, define \(\pmb {f}_\mathtt{t}(T)\), \(\mathtt{t}\in \{{\text {in}},\text {ex}\}\), to be the vector \(\pmb {w}\in {\mathbb {Z}}_+^{\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}\) that consists of the following entries:

• \(\pmb {w}({\mathtt{a}})=|\{v\in {\widetilde{V}}_\mathtt{t} \mid \alpha (v)={\mathtt{a}}\}|\), \({\mathtt{a}}\in \Lambda \),

• \(\pmb {w}(\gamma )=|\{uv\in {\widetilde{E}}_\mathtt{t} \mid \{\alpha (u),\alpha (v)\}=\{\mathtt{a,b}\}, \beta (uv)=q\}|\), \(\gamma =(\mathtt{a,b}, q)\in \Gamma \),

• \(\pmb {w}(\mu )=|\{uv\in {\widetilde{E}}_\mathtt{t}\mid \{\deg _T(u),\deg _T(v)\}=\{ d,d'\}, \beta (uv)=m\}|\), \(\mu =(d,d', m)\in \text {Bc}\),

• \(\pmb {w}(\text {dg}i)=|\{v\in {\widetilde{V}}_\mathtt{t}\mid \deg _T(v)=i\}|\), \(\text {dg}i\in \text {Dg}\).

Define \(\pmb {f}(T)\triangleq (\pmb {f}_{\text {in}}(T),\pmb {f}_\text {ex}(T))\). The entry for an element \(\mathtt{e}\in {\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}\) in \(\pmb {f}_\mathtt{t}(T)\), \(\mathtt{t}\in \{{\text {in}},\text {ex}\}\) is denoted by \(\pmb {f}_\mathtt{t}(\mathtt{e}; T)\). For a subset B of \( {\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}\), let \(\pmb {f}_{\mathtt{t}[B]}(T)\) denote the projection of \(\pmb {f}_\mathtt{t}(T)\) onto B.

Our aim is to generate all chemical bi-rooted (resp., tri-rooted) trees T with diameter \(\text {dia}^*\) such that \(\pmb {f}(T)=\pmb {x}^*\).

A new algorithm for computing chemical bi-rooted trees G with \({\text {bl}}_2(G)=2\)

This section describes a sketch of our new graph search algorithm for the case of \({\text {bl}}_2(G)=2\). See Appendix “A sketch of algorithm for computing chemical tri-rooted trees G with \({\text {bl}}_2(G)=3\)” for a sketch of a new algorithm for the case of \({\text {bl}}_2(G)=3\).

We call a chemical graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) with diameter \(\text {dia}^*\) and \({\text {bl}}_2(G)=2\) a target graph.

A chemical acyclic graph G with \({\text {bl}}_2(G)=2\) has exactly two leaf 2-branches \(v_i\), \(i=1,2\), where the length of the path between the two leaf 2-branches \(v_1\) and \(v_2\) of a target graph G is \(\text {dia}^*-2k^*=\text {dia}^*-4\). We observe that a connected subgraph T of a target graph G that satisfies (76) for \({\text {bl}}_2(G)=2\) is a chemical rooted or bi-rooted tree with roots u and v, where possibly \(u = v\). We call such a subgraph T an internal-subtree (resp., end-subtree) of G if neither (resp., one) of u and v is a 2-branch in G. When \(u=v\), we call an internal-subtree (resp., end-subtree) T of G an internal-fringe-tree (resp., end-fringe-tree) of  G. Figure 10a–d illustrates an internal-subtree, an internal-fringe-tree, an end-subtree and an end-fringe-tree of G.

An illustration of subtrees T of a chemical acyclic graph G in Fig. 6a, where the vertices/edges in T are depicted by solid lines: a An internal-subtree T of G; b An internal-fringe-tree T of G; c An end-subtree T of G; d An end-fringe-tree T of G

Full size image

Let \(\delta _1=\lfloor \frac{\text {dia}^* - 5}{2} \rfloor \) and \(\delta _2 =\text {dia}^*-5-\delta _1=\lceil \frac{\text {dia}^* - 5}{2} \rceil \). We regard a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) with \({\text {bl}}_2(G)=2\) and diameter \(\text {dia}^*\) as a combination of two chemical bi-rooted trees \(T_1\) and \(T_2\) with \(\ell (P_{T_i})=\delta _i\), \(i=1,2\), joined by an edge \(e=r_1(T_1)r_1(T_2)\), as illustrated in Fig. 11.

An illustration of combining two bi-rooted trees \(T_1=T_{\pmb {w}^1}\) and \(T_2=T_{\pmb {w}^2}\) with a new edge with multiplicity m joining vertices \(r_1(T_1)\) and \(r_1(T_2)\) to construct a target graph G, where \({\mathtt{a}}_i\in \Lambda \), \(d_i\in [1,d_\text {max}-1]\), \(m_i\in [d_i,{\text {val}}({\mathtt{a}}_i)-1]\), \(i = 1, 2\), and \(m \in [1, \min \{3, {\text {val}}({\mathtt{a}}_1) - m_1, {\text {val}}({\mathtt{a}}_2) - m_2\} ]\)

Full size image

We start with generating chemical rooted trees and then iteratively extend chemical bi-rooted trees T with \(\ell (P_T)=1,2,\dots ,\delta _1\), before we finally combine two chemical bi-rooted trees \(T_1\) and \(T_2\) with \(\ell (P_{T_i})=\delta _i\). To describe our algorithm, we introduce some notation.

• Let \({\mathcal {T}}(\pmb {x}^*)\) denote the set of all bi-rooted trees T (where possibly \(r_1(T)=r_2(T)\)) such that \(\pmb {f}_{{\text {in}}}(T)\le \pmb {x}_{{\text {in}}}^*\) and \(\pmb {f}_{\text {ex}}(T)\le \pmb {x}_{\text {ex}}^*\), which is a necessary condition for T to be an internal-subtree or end-subtree of a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\).

• Let \({\mathcal {FT}}\) denote the set of all rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) that can be a 2-fringe-tree of a target graph G, where T satisfies the size constraint (1) of 2-fringe-trees.

• For each integer \(h\in [1,\text {dia}^*-4]\), let \({\mathcal {T}}_{\text {end}}^{(h)}\) denote the set of all bi-rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) that can be an end-subtree of a target graph G such that \(\ell (P_T)=h\), and each 2-fringe-tree \(T_v\) rooted at a vertex v in \(P_T\) belongs to \({\mathcal {FT}}\).

The idea of our new algorithm is to compute only the set \(\text {W}_{\text {end}}^{(h)}\) of frequency vectors \(\pmb {w}\) of end trees, whose size \(|\text {W}_{\text {end}}^{(h)}|\) is much more restricted than that of \({\mathcal {T}}_{\text {end}}^{(h)}\). We compute the set \(\text {W}_{\text {end}}^{(h)}\) of frequency vectors \(\pmb {w}\) of trees in \({\mathcal {T}}_{\text {end}}^{(h)}\) iteratively for each integer \(h\ge 0\). During the computation, we keep a sample of a tree \(T_{\pmb {w}}\) for each frequency vector \(\pmb {w}\) so that a final step can construct some number of target graphs G by assembling these sample trees. Based on this, we generate target graphs \(G\in {\mathcal {G}}(\pmb {x}^*)\) by the following steps:

1. 1.
   1. (i)

      Compute \({\mathcal {FT}}\) by a branch-and-bound procedure that generates all possible rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) (where \(r_1(T)=r_2(T)\)) that can be a 2-fringe-tree of a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\);

   2. (ii)

      Compute the set \(\text {W}^{(0)}\) of all vectors \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\) for some tree \(T\in {\mathcal {FT}}\), and let \(\text {W}_{\text {end}}^{(0)} \subseteq \text {W}^{(0)}\) be those trees with height exactly 2;

   3. (iii)

      For each vector \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\in \text {W}^{(0)}\), choose a sample tree \(T_{\pmb {w}}\in {\mathcal {FT}}\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\), and store these sample trees;

2. 2.

   For each integer \(h=1,2,\ldots , \delta _2\), iteratively execute the next:

   1. (i)

      Compute the set \(\text {W}_{\text {end}}^{(h)}\) of all vectors \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\) for some bi-rooted tree \(T\in {\mathcal {T}}_{\text {end}}^{(h)}\), where such a vector \(\pmb {w}\) is obtained from a combination of vectors \(\pmb {w}'\in \text {W}^{(0)}\) and \(\pmb {w}''\in \text {W}_{\text {end}}^{(h-1)}\);

   2. (ii)

      For each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(h)}\), store a sample tree \(T_{\pmb {w}}\), which is obtained from a combination of sample trees \(T_{\pmb {w}'}\) with \(\pmb {w}'\in \text {W}^{(0)}\) and \(T_{\pmb {w}''}\) with \(\pmb {w}''\in \text {W}_{\text {end}}^{(h-1)}\);

3. 3.

   We call a pair of vectors \(\pmb {w}^1\in \text {W}_{\text {end}}^{(\delta _1)}\) and \(\pmb {w}^2\in \text {W}_{\text {end}}^{(\delta _2)}\) feasible, if it admits a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) such that \(\pmb {w}_{{\text {in}}}^1+\pmb {w}_{{\text {in}}}^2\le \pmb {x}^*_{{\text {in}}}\) and \(\pmb {w}_{\text {ex}}^1+\pmb {w}_{\text {ex}}^2\le \pmb {x}^*_{\text {ex}}\). Find the set \(\text {W}_\text {pair}\) of all feasible pairs of vectors \(\pmb {w}^1\) and \(\pmb {w}^2\);

4. 4.

   For each feasible vector pair \((\pmb {w}^1,\pmb {w}^2)\in \text {W}_\text {pair}\), construct a corresponding target graph G by combining the corresponding samples trees \(T_{\pmb {w}^1}\) and \(T_{\pmb {w}^2}\), as illustrated in Fig. 11.

Detailed descriptions of the five steps in the above algorithm can be found in Appendix “Case of two leaf 2-branches”.

For a relatively large instance with \(n^*\ge 40\) and \(\text {dia}^*\ge 20\), the number \(|\text {W}_\text {pair}|\) of feasible vector pairs in Step 4 is still very large. In fact, the size \(|\text {W}_{\text {end}}^{(h)}|\) of a vector set \(\text {W}_{\text {end}}^{(h)}\) to be computed in Step 2 can also be considerably large during an execution of the algorithm. For such a case, we impose a time limitation on the running time for computing \(\text {W}_{\text {end}}^{(h)}\) and a memory limitation on the number of vectors stored in a vector set \(\text {W}_{\text {end}}^{(h)}\). With these limitations, we can compute only a limited subset \({\widehat{\text {W}}}_{\text {end}}^{(h)}\) of each vector set \(\text {W}_{\text {end}}^{(h)}\) in Step 2. Even with such a subset \({\widehat{\text {W}}}_{\text {end}}^{(h)}\), we still can find a large size of a subset \({\widehat{\text {W}}}_\text {pair}\) of \(\text {W}_\text {pair}\) in Step 3.

Our algorithm also delivers a lower bound on the number of all target graphs \(G\in {\mathcal {G}}(\pmb {x}^*)\) in the following way. In Step 1, we also compute the number \(t(\pmb {w})\) of trees \(T\in {\mathcal {FT}}\) such that \(\pmb {w}=\pmb {f}(T)\) for each \(\pmb {w}\in \text {W}^{(0)}\). In Step 2, when a vector \(\pmb {w}\) is constructed from two vectors \(\pmb {w}'\) and \(\pmb {w}''\), we iteratively compute the number \(t(\pmb {w})\) of trees T such that \(\pmb {w}=\pmb {f}(T)\) by \(t(\pmb {w}):=t(\pmb {w}')\times t(\pmb {w}'')\). In Step 3, when a feasible vector pair \((\pmb {w}^1,\pmb {w}^2)\in \text {W}_\text {pair}\) is obtained, we know that the number of the corresponding target graphs G is \(t(\pmb {w}^1)\times t(\pmb {w}^2)\). Possibly we compute a subset \({\widehat{\text {W}}}_\text {pair}\) of \(\text {W}_\text {pair}\) in Step 3. Then \((1/2)\sum _{(\pmb {w}^1,\pmb {w}^2)\in {\widehat{\text {W}}}_\text {pair}}t(\pmb {w}^1)\times t(\pmb {w}^2)\) gives a lower bound on the number of target graphs \(G\in {\mathcal {G}}(\pmb {x}^*)\), where we divided by 2 since an axially symmetric target graph G can correspond to two vector pairs in \(\text {W}_\text {pair}\).

A sketch of algorithm for computing chemical tri-rooted trees G with \({\text {bl}}_2(G)=3\)

We call a chemical graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) with diameter \(\text {dia}^*\) and \({\text {bl}}_2(G)=3\) a target graph. Let \(n^*_{\text {inl}}\triangleq \sum _{{\mathtt{a}}\in \Lambda } \pmb {x}^*_{{\text {in}}}({\mathtt{a}})\), which is the number of 2-internal vertices in a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\).

A chemical acyclic graph G with \({\text {bl}}_2(G)=3\) has exactly three leaf 2-branches \(v_i\), \(i=1, 2, 3\), and exactly one 2-internal vertex \(v_4\) adjacent to three 2-internal vertices \(v'_i\), \(i=1, 2, 3\), as illustrated in Fig. 6(b). We call vertex \(v_4\) the joint-vertex of G. Without loss of generality assume that the length of the path \(P_{v_1,v_2}\) between \(v_1\) and \(v_2\) is \(\text {dia}^*-4\) and that the length of the path \(P_{v_1,v'_1}\) is not smaller than that of \(P_{v_2,v'_2}\).

Analogously with the case of \({\text {bl}}_2(G)=2\), we define internal-subtree (resp., end-subtree, internal-fringe-tree, and end-fringe-tree) of G, to be a connected subgraph \(G'\) that satisfies (76). Observe that G can be partitioned into three end-subtrees \(T_i\), \(i=1, 2, 3\), the 2-fringe-tree \(T_4\) rooted at the joint-vertex \(v_4\) and three edges \(v'_i v_4\), \(i=1,2,3\), where the backbone path \(P_{T_i}\) connects leaf 2-branch \(v_i\) and vertex \(v'_i\). In particular, we call the end-subtree of G that consists of \(T_1\), \(T_2\), \(T_4\), and edges \(v'_i v_4\), \(i=1,2\), the main-subtree of G, which consists of the path \(P_{v_1, v_2}\) and all the 2-fringe-trees rooted at vertices in \(P_{v_1, v_2}\). We call \(T_3\) the co-subtree of G.

Let \(\delta _i\), \(i=1,2,3\) denote the length of the backbone path of \(T_i\). Note that

• \(\delta _1+\delta _2+2=\text {dia}^*-4\) and \( \delta _1\ge \delta _2\ge \delta _3= n^*_{\text {inl}} - \text {dia}^* + 2\),

from which it follows that

• \(\delta _2 \in [\delta _3, \lfloor \text {dia}^*/2 \rfloor -3 ] \) and \(\delta _1 \in [\lceil \text {dia}^*/2\rceil -3, \text {dia}^* - 6 -\delta _3]\).

We regard a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) with \({\text {bl}}_2(G)=3\) and diameter \(\text {dia}^*\) as a combination of the main-subtree and the co-subtree joined with an edge. We represent the co-subtree as a chemical bi-rooted tree T with \(\ell (P_T)=\delta _3\). We represent the main-subtree of a target graph G as a tri-rooted tree T with \(\ell (P_T)=\text {dia}-4\) so that terminals \(r_1(T)\), \(r_2(T)\), and \(r_3(T)\), correspond to the two leaf 2-branches and the joint-vertex of G, respectively.

We start with generating chemical rooted trees and then iteratively extend chemical bi-rooted trees T with \(\ell (P_T)=1, 2, \dots , \text {dia}^*-6-\delta _3\), before we combine two chemical bi-rooted trees \(T'\) and \(T''\) to obtain a chemical tri-rooted tree \(T_{1}\) with \(\ell (P_{T_{1}})=\delta _i\), and finally, combine a chemical tri-rooted tree \(T_{1}\) and a chemical bi-rooted tree \(T_2\) with \(\ell (P_{T_2})=\delta _3\), to obtain a target graph \(G \in {\mathcal {G}}(\pmb {x}^*)\).

Analogously with the case of \({\text {bl}}_2(G)=2\), we define the set \({\mathcal {T}}(\pmb {x}^*)\) of all bi-rooted trees T, the set \({\mathcal {FT}}\) of all rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) that can be a 2-fringe-tree of a target graph G and the set \({\mathcal {T}}_{\text {end}}^{(h)}\), \(h\in [1, \text {dia}^* -6 -\delta _3]\), of all bi-rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) that can be an end-subtree of a target graph G such that \(\ell (P_T)=h\).

We generate target graphs \(G\in {\mathcal {G}}(\pmb {x}^*)\) by the following steps:

1. 1.

   Analogously with Step 1 for the case of \({\text {bl}}_2(G)=2\), compute the set \({\mathcal {FT}}\) by a branch-and-bound algorithm as described in "Step 1: Enumeration of 2-fringe-trees" section, and the set \(\text {W}^{(0)}\) of all vectors \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\) for some tree \(T\in {\mathcal {FT}}\). For each vector \(\pmb {w}\in \text {W}^{(0)}\), store a sample tree \(T_{\pmb {w}}\in {\mathcal {FT}}\), and let \(\text {W}_{\text {end}}^{(0)} \subseteq \text {W}^{(0)}\) be the set of feature vectors of possible end-trees with height 2;

2. 2.

   For each integer \(h=1,2,\ldots , \text {dia}^* -6 -\delta _3\), compute the set \(\text {W}_{\text {end}}^{(h)}\) of all vectors \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\) for some bi-rooted tree \(T\in {\mathcal {T}}_{\text {end}}^{(h)}\). For each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(h)}\), store a sample tree \(T_{\pmb {w}}\);

3. 3.

   For each integer \(h\in [\lceil \text {dia}^*/2\rceil -2, \text {dia}^* - 5 -\delta _3] \), compute the set \(\text {W}_{\text {end}+2}^{(h)}\) of all vectors \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\) of some bi-rooted tree T with \(\ell (P_T)=h\) that represents an end-subtree rooted at the joint-vertex. For each vector \(\pmb {w}\in \text {W}_{\text {end}+2}^{(h)}\), store a sample tree \(T_{\pmb {w}}\);

4. 4.

   For each integer \(\delta _1 \in [\lceil \text {dia}^*/2\rceil -3, \text {dia}^* - 6 -\delta _3]\), compute the set \(\text {W}_{\text {main}}^{(\delta _1+1)}\) of all vectors \(\pmb {w}=(\pmb {w}_{\text {in}},\pmb {w}_\text {ex})\) such that \(\pmb {w}_{\text {in}}=\pmb {f}_{{\text {in}}}(T)\) and \(\pmb {w}_\text {ex}=\pmb {f}_{\text {ex}}(T)\) for some tri-rooted tree T that represents the main-subtree such that the length of the path \(P_{r_2(T),r_3(T)}\) between terminals \(r_2(T)\) and \(r_3(T)\) is \(\delta _1+1\). For each vector \(\pmb {w}\in \text {W}_{\text {main}}^{(\delta _1+1)}\), store a sample tree \(T_{\pmb {w}}\);

5. 5.

   We call a pair of vectors \(\pmb {w}^1\in \text {W}_{\text {main}}^{(\delta _1+1)}\) and \(\pmb {w}^2\in \text {W}_{\text {end}}^{(\delta _3)}\) feasible if it admits a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) such that \(\pmb {w}_{{\text {in}}}^1+\pmb {w}_{{\text {in}}}^2\le \pmb {x}^*_{{\text {in}}}\) and \(\pmb {w}_{\text {ex}}^1+\pmb {w}_{\text {ex}}^2\le \pmb {x}^*_{\text {ex}}\). Find the set \(\text {W}_\text {pair}\) of all feasible pairs of vectors \(\pmb {w}^1\) and \(\pmb {w}^2\);

6. 6.

   For each feasible vector pair \((\pmb {w}^1,\pmb {w}^2)\in \text {W}_\text {pair}\), construct a corresponding target graph G by combining the samples trees \(T_{\pmb {w}^1}\) and \(T_{\pmb {w}^2}\), which correspond to the main-subtree and the co-subtree of a target graph G, respectively, as illustrated in Fig. 12.

An illustration of combining a tri-rooted \(T_1=T_{\pmb {w}^1}\) and a bi-rooted tree \(T_2=T_{\pmb {w}^2}\) with a new edge joining vertices \(r_3(T_1)\) and \(r_1(T_2)\) to construct a target graph G

Full size image

Detailed descriptions of the six steps in the above algorithm can be found in Appendix “Case of three leaf 2-branches”.

Frequency vectors of fictitious trees

Let T be a chemical bi-rooted or tri-rooted tree, where we regard a rooted tree T as a bi-rooted tree with \(r_1(T)=r_2(T)\) for a notational convenience. Recall that our algorithm generates a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) as a supergraph of T, where one of terminals \(r_1(T)\) and \(r_2(T)\) can be a 2-branch of G. We assume that the second terminal \(r_2(T)\) will be a 2-branch of G in such a case in our algorithms.

An illustration of fictitious rooted trees \(T[+p]\), \(p\in [1,3]\) for rooted trees T with \(r=r_1(T)=r_2(T)\) and \(d=\deg _T(r)\), where a dashed line depicts a fictitious edge incident to the terminal \(r_1(T)=r_2(T)\): (a) \(T[+1]\) and \(d=1\); (b) \(T[+1]\) and \(d=2\); (c) \(T[+1]\) and \(d=3\); (d) \(T[+2]\) and \(d=0\); (e) \(T[+2]\) and \(d=1\); (f) \(T[+2]\) and \(d=2\); (g) \(T[+3]\) and \(d=0\); (h) \(T[+3]\) and \(d=1\)

Full size image

An illustration of fictitious trees \(T[+q]\) and \(T\langle +1\rangle \) for bi-rooted trees and tri-rooted trees T: a \(T[+q]\) of a bi-rooted tree T; b \(T\langle +1\rangle \) of a tri-rooted tree T

Full size image

For an integer \(p\in [1,3]\), let \(T[+p]\) denote a fictitious chemical graph obtained from T by regarding the degree of terminal \(r_1(T)\) as \(\deg _T(r_1(T))+p\). Figure 13 (resp., Fig. 14a) illustrates fictitious trees \(T[+p]\) in the case of \(r_1(T)=r_2(T)\) (resp., \(r_1(T)\ne r_2(T)\)). The frequency vectors \(\pmb {f}_{{\text {in}}}(T[+p])\) and \(\pmb {f}_{\text {ex}}(T[+p])\) are obtained as follows: Let \(d=\deg _T(r_1(T))\), \(v_i\), \(i\in [1,d]\), denote the neighbors of \(r_1(T)\), and \(d_i=\deg _T(v_i)\), \(m_i=\beta (r_1(T) v_i)\), and \(\mu _i=(d, d_i, m_i)\), \(\mu '_i=(d+p, d_i, m_i)\), \(i\in [1,d]\).

For \(r_1(T)=r_2(T)\) and \(d'=d+p\),

• \(\pmb {f}_{{\text {in}}}(T[+p])= \pmb {f}_{{\text {in}}}(T) +\pmb {1}_{\text {dg}d'}-\pmb {1}_{\text {dg}d}\),   \(\displaystyle { \pmb {f}_{\text {ex}}(T[+p])= \pmb {f}_{\text {ex}}(T) + \sum _{1\le i\le d}( \pmb {1}_{\mu '_i} -\pmb {1}_{\mu _i}) }\).

For \(r_1(T)\ne r_2(T)\) and \(d'=d+p\), where \(v_d\) denotes the vertex in \(P_T\),

• \(\pmb {f}_{{\text {in}}}(T[+1]) = \pmb {f}_{{\text {in}}}(T) +\pmb {1}_{\text {dg}d'}-\pmb {1}_{\text {dg}d} +\pmb {1}_{\mu '_d}-\pmb {1}_{\mu _d}\),

• \(\displaystyle { \pmb {f}_{\text {ex}}(T[+1]) = \pmb {f}_{\text {ex}}(T) +\sum _{1\le i\le d-1}(\pmb {1}_{\mu '_i}-\pmb {1}_{\mu _i}) }\).

Let T be a chemical tri-rooted tree, where the third terminal \(r_3(T)\) is in the backbone path \(P_T\) between vertices \(r_1(T)\) and \(r_2(T)\). Let \(T\langle +1\rangle \) denote a fictitious chemical graph obtained from T by regarding the degree of terminal \(r_3(T)\) as \(\deg _T(r_3(T))+1\). Figure 14b illustrates a fictitious tri-rooted tree \(T\langle +1\rangle \). The frequency vectors \(\pmb {f}_{{\text {in}}}(T\langle +1\rangle )\) and \(\pmb {f}_{\text {ex}}(T\langle +1\rangle )\) are obtained as follows: Let \(d=\deg _T(r_3(T))\), \(v_i\), \(i\in [1,d]\), denote the neighbors of \(r_3(T)\), where \(v_{d-1}\) and \(v_{d}\) are contained in the path \(P_T\). For each index \(i\in [1,d]\), let \(d_i=\deg _T(v_i)\), \(m_i=\beta (r_3(T) v_i)\), \(\mu _i=(d, d_i, m_i)\), and \(\mu '_i=(d+1, d_i, m_i)\).

Then

Sets of frequency vectors

For an element \({\mathtt{a}}\in \Lambda \) and integers \(d \in [0, d_\text {max}-2]\) and \(m\in [d,{\text {val}}({\mathtt{a}})-1]\), let \(\text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d, m)\) (resp., \(\text {W}_{\text {inl}+3}^{(0)}({\mathtt{a}}, d, m)\)) denote the set of frequency vectors \((\pmb {f}_{{\text {in}}}(T[+2]), \pmb {f}_{\text {ex}}(T[+2]))\) (resp., \((\pmb {f}_{{\text {in}}}(T[+3]), \pmb {f}_{\text {ex}}(T[+3]))\)) of a chemical rooted tree T such that

• \(r_1(T)=r_2(T)\), the height of T is at most 2,

• \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\), and \(\beta (r_1(T))=m\).

Recall that \(\beta (u)=\sum _{uv\in E}\beta (uv)\), defined in “Preliminary” section.

For an element \({\mathtt{a}}\in \Lambda \) and integers \(d \in [1, d_\text {max}-1]\), \(m\in [d,{\text {val}}({\mathtt{a}})-1]\), and \(h\ge 0\), let \(\text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d, m)\) (resp., \(\text {W}_{\text {end}+2}^{(h)}({\mathtt{a}}, d, m)\)) denote the set of frequency vectors \((\pmb {f}_{{\text {in}}}(T[+1]), \pmb {f}_{\text {ex}}(T[+1]))\) (resp., \((\pmb {f}_{{\text {in}}}(T[+2]), \pmb {f}_{\text {ex}}(T[+2]))\)) of chemical bi-rooted trees T such that

• \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\), \(\beta (r_1(T))=m\), \(\ell (P_T)=h\) and

• if \(h=0\) then the height of the tree \(T'\) rooted at \(r_2(T)\) is 2.

Case of two leaf 2-branches

Step 1: Enumeration of 2-fringe-trees

The main task of Step 1 is to compute for each tuple \(({\mathtt{a}},d,m)\) of an element \({\mathtt{a}}\in \Lambda \) and integers \(d \in [1,d_\text {max}-1]\) (resp., \(d\in [0,d_\text {max}-2]\)) and \(m\in [d,{\text {val}}({\mathtt{a}})-1]\) (resp., \(m\in [d,{\text {val}}({\mathtt{a}})-2]\)), the set \(\text {W}_{\text {end}}^{(0)}({\mathtt{a}}, d, m)\) (resp., \(\text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d, m)\)) of all frequency vectors \(\pmb {f}(T[+1])\) (resp., \(\pmb {f}(T[+2])\)) of chemical rooted trees T such that \(r_1(T)=r_2(T)\), \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\) and \(\beta (r_1(T))=m\).

Step 1 first computes the set \({\mathcal {FT}}\) of all possible chemical rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) (where \(r_1(T)=r_2(T)\)) that can be a 2-fringe-tree of a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\). For this, we design a branch-and-bound procedure where we append a new vertex one by one to construct a rooted tree with only one child. To design a bounding procedure, we derive a property of the structure of chemical rooted trees that can be a 2-fringe-tree of a target graph.

Let \(G_0\) be a chemical rooted tree with a terminal \(r_0=r_1(G_0)=r_2(G_0)\), where \(\pmb {f}_{{\text {in}}}(\alpha (r_0);G_0)=1\) and \(\pmb {f}_{{\text {in}}}({\mathtt{a}};G_0)=0\), \({\mathtt{a}}\in \Lambda \setminus \{\alpha (r_0)\}\) and \(\pmb {f}_{{\text {in}}}(\gamma ;G_0)=0\), \(\gamma \in \Gamma \). For a vector \(\pmb {x}= (\pmb {x}_{\text {in}}, \pmb {x}_\text {ex})\) with \(\pmb {x}_{\text {in}}, \pmb {x}_\text {ex}\in {\mathbb {Z}}_+^{ {\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}}\), we call \(G_0\) \(\pmb {x}\)-extensible if some chemical acyclic graph \(G\in {\mathcal {G}}(\pmb {x})\) contains \(G_0\) as a subgraph of a 2-fringe-tree T rooted at \(r_0\) in G.

We use the next condition as a bounding procedure when we generate chemical rooted trees in Step 1.

Lemma 3

For a branch-parameter \(k=2\), let \(\pmb {x}^* = (\pmb {x}_{\text {in}}^*, \pmb {x}_\text {ex}^*)\) be a vector with \(\pmb {x}_{\text {in}}^*, \pmb {x}_\text {ex}^*\in {\mathbb {Z}}_+^{{\Lambda \cup \Gamma \cup \text {Bc}\cup \text {Dg}}}\), and \(G_0\) be a chemical rooted tree rooted at a vertex \(r_0\) such that \(\pmb {f}(G_0)\le \pmb {x}^*\).

1. (i)

   Graph \(G_0\) is \(\pmb {x}^*\)-extensible only when the next holds for any subset \(\Lambda '\subseteq \Lambda \):

   $$\begin{aligned} \sum _{{\mathtt{a}}\in \Lambda '} (\pmb {x}_\text {ex}^*({\mathtt{a}}) - \pmb {f}_\text {ex}({\mathtt{a}}; G_0))&\le \sum _{\begin{array}{c} \gamma =(\mathtt{a,b},m)\in \Gamma : \\ {\mathtt{a}}\in \Lambda ', \mathtt{b}\in \Lambda \setminus \Lambda ' \end{array}} (\pmb {x}_\text {ex}^*(\gamma ) - \pmb {f}_\text {ex}(\gamma ; G_0)) \nonumber \\&\quad + 2 \sum _{\begin{array}{c} \gamma =(\mathtt{a,b},m)\in \Gamma : \\ \mathtt{a,b}\in \Lambda ' \end{array}} (\pmb {x}_\text {ex}^*(\gamma ) - \pmb {f}_\text {ex}(\gamma ; G_0)). \end{aligned}$$

   (78)
2. (ii)

   Let \(G_1\) denote the chemical rooted tree obtained from \(G_0\) by appending a new atom with an element \(\mathtt{b}\in \Lambda \) to an atom with an element \({\mathtt{a}}\in \Lambda \) in \(G_0\) with a multiplicity q; i.e., we join an atom \({\mathtt{a}}\) in \(G_0\) and a new atom \(\mathtt{b}\) with an adjacency-configuration \((\mathtt{a,b},q)\). Then \(G_1\) is \(\pmb {x}^*\)-extensible only when the next holds:

   • \(\pmb {x}_\text {ex}^*({\mathtt{a}}) - \pmb {f}_\text {ex}({\mathtt{a}}; G_0)\le \pmb {nb}({\mathtt{a}}) -1\)

   for

   • \(\pmb {nb}({\mathtt{a}})= \displaystyle { \sum _{\begin{array}{c} \gamma =(\mathtt{a,b},m)\in \Gamma : \\ \mathtt{b\ne a}\in \Lambda \end{array} } (\pmb {x}_\text {ex}^*(\gamma ) - \pmb {f}_\text {ex}(\gamma ; G_0) ) + 2 \sum _{\gamma =(\mathtt{a,a},m)\in \Gamma } (\pmb {x}_\text {ex}^*(\gamma ) - \pmb {f}_\text {ex}(\gamma ; G_0) ) }\).

Proof

1. (i)

   Assume that \(G_0\) is a subgraph of a 2-fringe-tree T in some chemical graph \(G\in {\mathcal {G}}(\pmb {x}^*)\) so that T is rooted at \(r_0\). The left-hand side means the number of the remaining 2-external vertices with elements in \(\Lambda '\) in the 2-fringe-trees in G. Each of such atoms has a neighbor in the connected graph G. The right-hand side indicates an upper bound on the number of 2-external edges joining elements in \(\Lambda '\) in the 2-fringe-trees in G.

2. (ii)

   Note that \(\pmb {f}_{\text {ex}[\Lambda \cup \Gamma ]}(G_1) = \pmb {f}_{\text {ex}[\Lambda \cup \Gamma ]}(G_0)+\pmb {1}_\mathtt{b}+\pmb {1}_{\gamma }\). For \(\Lambda '=\{{\mathtt{a}}\}\), the left-hand side in Eq. (78) is \(\pmb {x}^*_\text {ex}({\mathtt{a}})- \pmb {f}_\text {ex}({\mathtt{a}}; G_0)\), which remains unchanged if \(\mathtt{a\ne b}\) (resp., is reduced by 1 if \(\mathtt{a= b}\)); and the right-hand side in (78) is \(\pmb {nb}({\mathtt{a}})\), which is reduced by 1 if \(\mathtt{a\ne b}\) (resp., is reduced by 2 if \(\mathtt{a= b}\)). That is, the left-hand side minus the right-hand side in (78) is always reduced by 1. This gives the required necessary condition for \(G_1\) to be \(\pmb {x}^*\)-extensible.

\(\square \)

An illustration of rooted trees T with height at most 2 and only one child satisfying the size constraint: a case of \(n(T)=2\); b case of \(n(T)=3\); c case of \(n(T)=4\); d case of \(n(T)=5\)

Full size image

Figure 15 illustrates all graph structures of rooted trees T with height at most 2 and only one child satisfying the size constraint (1). For each element \({\mathtt{a}}\in \Lambda \), we enumerate chemical trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) rooted at vertex r with \(\alpha (r)={\mathtt{a}}\) that has only one child by a branch-and-bound algorithm. Let \({\mathcal {T}}_{\mathtt{a}}\) denote the set of resulting rooted trees for each root element \({\mathtt{a}}\in \Lambda \).

We next enumerate chemical trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) rooted at vertex r with \(\alpha (r)={\mathtt{a}}\) that has two or three children by generating a combination of two or three graphs in \({\mathcal {T}}_{\mathtt{a}}\). During generating graphs, our bounding procedure tests whether the current graph satisfies the necessary condition in Lemma 3(ii).

Finally, we compute the following sets:

for each element \({\mathtt{a}}\in \Lambda \), integers \(d\in [1,d_\text {max}-1]\), \(m\in [d,{\text {val}}({\mathtt{a}})-1]\), the set \(\text {W}_{\text {end}}^{(0)}({\mathtt{a}}, d, m)\) of frequency vectors \(\pmb {f}(T[+1])\) for rooted trees \(T\in {\mathcal {T}}_{\mathtt{a}}\) with \(\deg _T(r)=d\) and height 2;

for each element \({\mathtt{a}}\in \Lambda \), integers \(d\in [0,d_\text {max}-2]\), \(m\in [d,{\text {val}}({\mathtt{a}})-2]\), the set \(\text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d, m)\) of frequency vectors \(\pmb {f}(T[+2])\) for rooted trees \(T\in {\mathcal {T}}_{\mathtt{a}}\) with \(\deg _T(r)=d\) and height at most 2.

For each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(0)}({\mathtt{a}}, d, m)\) (resp., \(\pmb {w}\in \text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d, m)\)), we store a sample tree \(T_{\pmb {w}}\).

We remark that the size of the set \({\mathcal {FT}}\) depends on the vector \(\pmb {x}^*\). However, since the height of trees is limited to 2, the degree is at most 3 or 4, and the size constraint (1) on fringe trees in "Our target graph class" section, the size of the set \({\mathcal {FT}}\) is fairly limited.

Step 2: Generation of frequency vectors of end-subtrees

The main task of Step 2 is to compute the following sets in the ascending order of \(h=1,2,\ldots , \delta _2 \):

For elements \({\mathtt{a}}\in \Lambda \), integers \(d \in [1,d_\text {max}-1]\), \(m \in [d ,{\text {val}}({\mathtt{a}})-1]\), and \(h\in [ 1, \delta _2 ]\), the sets \(\text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d , m )\) of all frequency vectors \(\pmb {f}(T[+1])\) of chemical bi-rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) such that \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d \), \(\beta (r_1(T))=m\) and \(\ell (P_T)=h\).

Observe that each vector \(\pmb {w}=(\pmb {w}_{{\text {in}}},\pmb {w}_{\text {ex}})\in \text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d , m )\) is obtained from a combination of vectors \(\pmb {w}'=(\pmb {w}'_{{\text {in}}},\pmb {w}'_{\text {ex}})\in \text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d -1, m')\) and \(\pmb {w}''=(\pmb {w}''_{{\text {in}}},\pmb {w}''_{\text {ex}}) \in \text {W}_{\text {end}}^{(h-1)}(\mathtt{b}, d'', m'')\) such that

• \(m' \le {\text {val}}({\mathtt{a}}) - 2\), \(1\le m - m' \le {\text {val}}(b) -m'' \),

• \(\pmb {w}_{{\text {in}}} = \pmb {w}'_{{\text {in}}} + \pmb {w}''_{{\text {in}}} + \pmb {1}_\gamma + \pmb {1}_\mu \le \pmb {x}_{{\text {in}}}^*\), \(\pmb {w}_{\text {ex}} = \pmb {w}'_{\text {ex}} + \pmb {w}''_{\text {ex}} \le \pmb {x}_{\text {ex}}^*\)

• for \(\gamma = ({\mathtt{a}}, \mathtt{b}, m - m' ) \in \Gamma \) and \(\mu = (d+1, d''+1, m - m') \in \text {Bc}\).

Figure 16 illustrates this process of computing a vector \(\pmb {w}\in \text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d , m )\).

An illustration of appending a rooted tree \(T'\) to a bi-rooted tree \(T''\) to compute a vector \(\pmb {w}\in \text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d , m )\) from the frequency vectors \(\pmb {w}'=\pmb {f}(T'[+2])\in \text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d -1, m')\) of a rooted tree \(T'\) and \(\pmb {w}'' =\pmb {f}(T''[+1]) \in \text {W}_{\text {end}}^{(h-1)}(\mathtt{b}, d'', m'')\) of a bi-rooted tree \(T''\)

Full size image

For each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d , m )\) obtained from a combination \(\pmb {w}'\in \text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d -1, m')\) and \(\pmb {w}'' \in \text {W}_{\text {end}}^{(h-1)}(\mathtt{b}, d'', m'')\), we construct a sample tree \(T_{\pmb {w}}\) from their sample trees \(T_{\pmb {w}'}\) and \(T_{\pmb {w}''}\).

Step 3: Enumeration of feasible vector pairs

A feasible pair of vectors is defined to be a pair of vectors \(\pmb {w}^i= (\pmb {w}^i_{{\text {in}}}, \pmb {w}^i_{\text {ex}}) \in \text {W}_{\text {end}}^{(\delta _i)} ({\mathtt{a}}_i,d_i,m_i)\), \({\mathtt{a}}_i\in \Lambda \), \(d_i\in [1,d_\text {max}-1]\), \(m_i\in [d_i,{\text {val}}({\mathtt{a}}_i)-1]\), \(i = 1,2\) that admits an adjacency-configuration \(\gamma = ({\mathtt{a}}_1, {\mathtt{a}}_2, m) \in \Gamma \) and a bond-configuration \(\mu = (d_1 + 1, d_2 + 1, m) \in \text {Bc}\) with an integer \(m \in [1, \min \{3, {\text {val}}({\mathtt{a}}_1) - m_1, {\text {val}}({\mathtt{a}}_2) - m_2\} ]\) such that

• \(\pmb {x}^*_{{\text {in}}} = \pmb {w}^1_{{\text {in}}} + \pmb {w}^2_{{\text {in}}} + \pmb {1}_{\gamma } + \pmb {1}_{\mu }\) and \(\pmb {x}^*_{\text {ex}} = \pmb {w}^1_{\text {ex}} + \pmb {w}^2_{\text {ex}}\),

or equivalently \(\pmb {w}^1\) is equal to the vector \((\pmb {x}^*_{{\text {in}}} - \pmb {w}^2_{{\text {in}}}-\pmb {1}_{\gamma } -\pmb {1}_{\mu }, \pmb {x}^*_{\text {ex}} - \pmb {w}^1_{\text {ex}})\), which we call the \((\gamma , \mu )\)-complement of \(\pmb {w}^2\), and denote it by \(\overline{\pmb {w}^2}\).

The main task of Step 3 is to enumerate all feasible vector pairs \((\pmb {w}^1,\pmb {w}^2)\), \(\pmb {w}^i \in \text {W}_{\text {end}}^{(\delta _i)}({\mathtt{a}}_i,d_i,m_i)\) with \({\mathtt{a}}_i\in \Lambda \), \(d_i\in [1,d_\text {max}-1]\), \(m_i\in [d_i,{\text {val}}({\mathtt{a}}_i)-1]\), \(i = 1,2\).

To efficiently search for a feasible pair of vectors in two sets \(\text {W}_{\text {end}}^{(\delta _i)}({\mathtt{a}}_i,d_i,m_i)\), \(i = 1,2\), we first compute the \((\gamma , \mu )\)-complement vector \(\overline{\pmb {w}}\) of each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(\delta _2)}({\mathtt{a}}_2,d_2,m_2)\) for each pair of \( \gamma = ({\mathtt{a}}_1, {\mathtt{a}}_2, m) \in \Gamma \) and \(\mu = (d_1 + 1, d_2 + 1, m) \in \text {Bc}\) with \(m \in [1, \min \{3, {\text {val}}({\mathtt{a}}_1) - m_1, {\text {val}}({\mathtt{a}}_2) - m_2\} ]\), and denote by \(\overline{\text {W}_{\text {end}}^{(\delta _2)}}\) the set of the resulting \((\gamma , \mu )\)-complement vectors. Observe that \((\pmb {w}^1,\pmb {w}^2)\) is a feasible vector pair if and only if \(\pmb {w}_1=\overline{\pmb {w}_2}\). To find such pairs, we merge the sets \(\text {W}_{\text {end}}^{(\delta _1)}({\mathtt{a}}_1,d_1,m_1)\) and \(\overline{\text {W}_{\text {end}}^{(\delta _2)}}\) into a sorted list \(L_{\gamma , \mu }\). Then each feasible vector pair \((\pmb {w}^1,\pmb {w}^2)\) appears as a consecutive pair of vectors \(\pmb {w}_1\) and \(\overline{\pmb {w}_2}\) in the list \(L_{\gamma , \mu }\).

Step 4: Construction of chemical graphs

The task of Step 4 is to construct for each feasible vector pair \(\pmb {w}^i\in \text {W}_{\text {end}}^{(\delta _i)}({\mathtt{a}}_i,d_i,m_i)\), \(i = 1,2\) such that \(\pmb {w}^1\) is equal to the \((\gamma = ({\mathtt{a}}_1, {\mathtt{a}}_2, m), \mu )\)-complement vector \(\overline{\pmb {w}^2}\) of \(\pmb {w}^2\), construct a target graph \(T_{(\pmb {w}_1, \pmb {w}_2)}\in {\mathcal {G}}(\pmb {x}^*)\) by combining the sample trees \(T_i=T_{\pmb {w}^i}\) of vectors \(\pmb {w}^i\) with an edge \(e=r_1(T_1)r_1(T_2)\) such that \(\beta (e)=m\). Figure 11 illustrates two sample trees \(T_i\), \(i=1,2\) to be combined with a new edge \(e=r_1(T_1)r_1(T_2)\).

Case of three leaf 2-branches

Step 1: Enumeration of 2-fringe-trees

The main task of Step 1 is to compute the following sets:

for each tuple \(({\mathtt{a}},d,m)\) of an element \({\mathtt{a}}\in \Lambda \) and integers \(d \in [1,d_\text {max}-1]\) (resp., \(d\in [0,d_\text {max}-2]\) and \(d\in [0,d_\text {max}-3]\)) and \(m\in [d,{\text {val}}({\mathtt{a}})-1]\) (resp., \(m\in [d,{\text {val}}({\mathtt{a}})-2]\) and \(m\in [d,{\text {val}}({\mathtt{a}})-3]\)), the set \(\text {W}_{\text {end}}^{(0)}({\mathtt{a}}, d, m)\) (resp., \(\text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d, m)\) and \(\text {W}_{\text {inl}+3}^{(0)}({\mathtt{a}}, d, m)\)) of all frequency vectors \(\pmb {f}(T[+1])\) (resp., \(\pmb {f}(T[+2])\) and \(\pmb {f}(T[+3])\)) of chemical rooted trees T such that \(r_1(T)=r_2(T)\), \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\) and \(\beta (r_1(T))=m\). For each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(0)}({\mathtt{a}}, d, m)\) (resp., \(\pmb {w}\in \text {W}_{\text {inl}}^{(0)}({\mathtt{a}}, d, m)\) and \(\pmb {w}\in \text {W}_{\text {inl}+3}^{(0)}({\mathtt{a}}, d, m)\)), we store a sample tree \(T_{\pmb {w}}\). This step can be designed in a similar way as Step 1 for the case of \({\text {bl}}_2(G)=2\).

Step 2: Generation of frequency vectors of end-subtrees

Analogously with Step 2 for the case of \({\text {bl}}_2(G)=2\), Step 2 computes the following sets in the ascending order of \(h=1,2,\ldots , \text {dia}^* -6 -\delta _3\):

For elements \({\mathtt{a}}\in \Lambda \), integers \(d \in [1,d_\text {max}-1]\), \(m \in [d ,{\text {val}}({\mathtt{a}})-1]\), \(i=1,2\), and \(h\in [ 1, \text {dia}^* -6 -\delta _3]\), the sets \(\text {W}_{\text {end}}^{(h)}({\mathtt{a}},d,m)\) of all frequency vectors \(\pmb {f}(T[+1])\) of chemical bi-rooted trees \(T\in {\mathcal {T}}(\pmb {x}^*)\) such that \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\), \(\beta (r_1(T))=m\) and \(\ell (P_T)=h\).

For each vector \(\pmb {w}\in \text {W}_{\text {end}}^{(h)}({\mathtt{a}}, d, m)\), we construct a sample tree \(T_{\pmb {w}}\) from their sample trees \(T_{\pmb {w}'}\) and \(T_{\pmb {w}''}\).

Step 3: Generation of frequency vectors of end-subtrees with two fictitious edges

The main task of Step 3 is to compute the following sets:

For elements \({\mathtt{a}}\in \Lambda \), integers \(d \in [1,d_\text {max}-2]\), \(m\in [d,{\text {val}}({\mathtt{a}})-2]\) and \(h\in [\lceil \text {dia}^*/2\rceil -2, \text {dia}^* - 5 -\delta _3] \), the sets \(\text {W}_{\text {end}+2}^{(h)}({\mathtt{a}},d,m)\) of all frequency vectors of bi-rooted trees \(T[+2]\) such that \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\), \(\beta (r_1(T))=m\) and \(\ell (P_T)=h\). For each vector \(\pmb {w}\in \text {W}_{\text {end}+2}^{(h)}({\mathtt{a}},d,m)\), we store a sample tree \(T_{\pmb {w}}\). This step can be designed in a similar way as Step 3 for the case of \({\text {bl}}_2(G)=2\).

Step 4: Enumeration of frequency vectors of main-subtrees

For an element \({\mathtt{a}}\in \Lambda \), and integers \(d \in [2,d_\text {max}-1]\), \(m \in [d, {\text {val}}({\mathtt{a}}) -1]\), and \(\delta _1 \in [\lceil \text {dia}^*/2\rceil -3, \text {dia}^* - 6 -\delta _3]\), define \(\text {W}_{\text {main}}^{(\delta _1+1)}({\mathtt{a}}, d, m)\) to be the set of the frequency vectors \(\pmb {f}(T\langle +1\rangle )\) of chemical tri-rooted trees T such that

• \(\alpha (r_1(T))={\mathtt{a}}\), \(\deg _T(r_1(T))=d\), \(\beta (r_1(T))=m\), \(\ell (P_T)=\text {dia}^*-4\) and

• the length of the path \(P_{r_2(T),r_3(T)}\) between vertices \(r_2(T)\) and \(r_3(T)\) is \(\delta _1+1\).

See Fig. 12 for the structure of a main-tree. Such a chemical tri-rooted graph T corresponds to the main-subtree of a target graph \(G\in {\mathcal {G}}(\pmb {x}^*)\).

The main task of Step 4 is to compute the sets \(\text {W}_{\text {main}}^{(\delta _1+1)}({\mathtt{a}}, d, m)\), \({\mathtt{a}}\in \Lambda \), \(d \in [2,d_\text {max}-1]\), \(m \in [d, {\text {val}}({\mathtt{a}}) -1]\), \(\delta _1 \in [\lceil \text {dia}^*/2\rceil -3, \text {dia}^* - 6 -\delta _3]\). Each vector \(\pmb {w}\in \text {W}_{\text {main}}^{(\delta _1+1)}({\mathtt{a}}, d, m)\) can be obtained from a combination of vectors \(\pmb {w}^1\in \text {W}_{\text {end}+2}^{(\delta _1+1)}({\mathtt{a}}, d-1, m'')\) and \(\pmb {w}^2\in \text {W}_{\text {end}}^{(\delta _2)}({\mathtt{a}}', d', m')\) such that \(\delta _1 + \delta _2 = \text {dia}^* - 4\) and \(\delta _1 \ge \delta _2 \), as illustrated in Fig. 17. For each vector \(\pmb {w}\in \text {W}_{\text {main}}^{(\delta _1+1)}({\mathtt{a}}, d, m)\), we store a sample tree \(T_{\pmb {w}}\). This step can be designed in a similar way as Step 3 for the case of \({\text {bl}}_2(G)=2\).

An illustration of computing the frequency vector \(\pmb {w}=\pmb {f}(T\langle +1\rangle ) \in \text {W}_{\text {main}}^{(\delta _1+1)}({\mathtt{a}}, d, m)\) of a tri-rooted tree T from the frequency vectors \(\pmb {w}^1=\pmb {f}(T_1[+2])\in \text {W}_{\text {end}+2}^{(\delta _1+1)}({\mathtt{a}}, d-1, m'')\) and \(\pmb {w}^2=\pmb {f}(T_2[+1])\in \text {W}_{\text {end}}^{(\delta _2)}({\mathtt{a}}', d', m')\) for bi-rooted trees \(T_1\) and \(T_2\)

Full size image

Step 5: Enumeration of feasible vector pairs

Analogously with the case of \({\text {bl}}_2(G)=2\), a feasible pair of vectors is defined to be a pair of vectors \(\pmb {w}^1= (\pmb {w}^1_{{\text {in}}}, \pmb {w}^1_{\text {ex}}) \in \text {W}_{\text {main}}^{(\delta _1+1)}({\mathtt{a}}_1,d_1, m_1)\), and \(\pmb {w}^2= (\pmb {w}^2_{{\text {in}}}, \pmb {w}^2_{\text {ex}}) \in \text {W}_{\text {end}}^{(\delta _3)}({\mathtt{a}}_2,d_2, m_2)\), \(\delta _1 \in [\lceil \text {dia}^*/2\rceil -3, \text {dia}^* - 6 -\delta _3]\), \({\mathtt{a}}_i\in \Lambda \), \(d_i\in [1,d_\text {max}-1]\), \(m_i\in [d_i,{\text {val}}({\mathtt{a}}_i)-1]\), \(i = 1,2\) that admits an adjacency-configuration \(\gamma = ({\mathtt{a}}_1, {\mathtt{a}}_2, m) \in \Gamma \) and a bond-configuration \(\mu = (d_1 + 1, d_2 + 1, m) \in \text {Bc}\) with an integer \(m \in [1, \min \{3, {\text {val}}({\mathtt{a}}_1) - m_1, {\text {val}}({\mathtt{a}}_2) - m_2\} ]\) such that

• \(\pmb {x}^*_{{\text {in}}} = \pmb {w}^1_{{\text {in}}} + \pmb {w}^2_{{\text {in}}} + \pmb {1}_{\gamma } + \pmb {1}_{\mu }\) and \(\pmb {x}^*_{\text {ex}} = \pmb {w}^1_{\text {ex}} + \pmb {w}^2_{\text {ex}}\).

Step 5 computes the set of all feasible vector pairs \((\pmb {w}^1,\pmb {w}^2)\) by using a sorting algorithm as in the Step 4 for the case of \({\text {bl}}_2(G)=2\).

Step 6: Construction of chemical graphs

Analogously with Step 4 for the case of \({\text {bl}}_2(G)=2\), Step 6 constructs a target graph \(T_{(\pmb {w}_1, \pmb {w}_2)}\in {\mathcal {G}}(\pmb {x}^*)\) for each feasible vector pair \((\pmb {w}^1,\pmb {w}^2)\) by combining the sample trees \(T_i=T_{\pmb {w}^i}\) of vectors \(\pmb {w}^i\) with a new edge \(e=r_1(T_1)r_1(T_2)\).

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions