Vyxal, 14 bytes

Þ∪Þx:ƛ□꘍GN;ÞMİ

Try it Online!

I think Þ∪ can be J, but that makes it so slow that it's practically impossible to run.

This is already really, really, really slow.

Þ∪Þx:ƛ□꘍GN;ÞMİ
Þ∪             # Get the multiset union of the two strings
  Þx           # Get all combinations without replacement of all lengths and of all orders
    :          # Duplicate (call this X)
     ƛ         # Over each:
      □꘍       #  For each of the inputs, get the Levenshthein distance between that input and this string
        G      #  Get the maximum of that
         N     #  Negate that
          ;    # Close map
           ÞM  # Get the indices of maximal elements. (If Vyxal had an "indices of minimal elements" builtin, then negate could be removed.)
             İ # Index this into X

05AB1E, 9 bytes

JæΣ.L{R}н

Outputs one possible result.

Try it online or verify all test cases.

Explanation:

J        # Join the (implicit) input-pair together to a single string
 æ       # Pop and get the powerset of this string
  Σ      # Sort it by:
   .L    #  Get the Levenshtein distance of the two strings with the (implicit) input
     {   #  Sort the pair from lowest to highest
      R  #  Reverse it from highest to lowest
  }н     # After the sort-by: leave the first result
         # (which is output implicitly as result)

Note: if we'd use à (max) instead of {R, it would also give results where the Levensthein distance is [2,2] even if [1,2] are available. If we'd use just {, it would give results where the Levensthein distance is [0,3] over [1,2]. So using the {R ensures the minimum pair while prioritizing the maximum of the pair.

Although it works for the test cases, I'm not 100% sure if Jæ is always correct. Using an additional group-by (and uniquify) to get all results shows that Stegfan and Segfan aren't among them for example, since the powerset builtin uses characters of the joined string in order. This could be fixed by changing it to Jœ€æ˜ [+3 bytes] (powerset of each permutation of the joined string), but it becomes extremely slow and it'll timeout.

Factor + lcs math.combinatorics, ^{123 120} 78 bytes

[ 2dup append all-subsets [ tuck [ levenshtein ] 2bi@ max ] 2with infimum-by ]

Attempt This Online!

-42 bytes(!) using the idea behind @KevinCruijssen's 05AB1E answer.

Vyxal, 9 bytes

∑ṗµ?꘍sṘ;h

Try it Online!

Port of Kevin Cruijssen's answer, go upvote that!

Pyth, 95 93 bytes

J._+k._KwFGUz=ZJ=J]+_._<zhGkVlK aJh.mlb[+<zhG@ZhN+eJ<KhN+m+d@zG@ZN?q@zG@KNY]+e@ZN@KN;@eJ/leJ2

Try it online!

I know what you're thinking: "That seems long." In my defense, it's been a while since I've touched Pyth. Furthermore, this runs in Ο(m·n·max(m,n)) time, where m, n are the lengths of the two input strings. This is wicked fast and can handle pretty long strings with ease. I used this algorithm from Wikipedia, but instead of keeping track of the distance, I keep track of the optimal path throughout. Then at the end we can simply output the middle element of the optimal path. This must satisfy the requirements, as if there were some string with a lower max distance, there would then exist a shorter path.

Python, 150 bytes

n=lambda A,B:{a[:x]+y+a[z:]for a in A for x in range(len(a)+1)for z in[x,x+1]for b in B for y in[*b,""]}
f=lambda A,B:(N:=n(A,B))&(M:=n(B,A))or f(N,M)

Attempt This Online!

Wolfram Language (Mathematica), 92 88 bytes

(c=Characters/@#;""<>#&/@MinimalBy[Subsets@*Join@@c,Max[EditDistance[x,#]&/@c]/.x->#&])&

Try it online!